home *** CD-ROM | disk | FTP | other *** search
/ Aminet 30 / Aminet 30 (1999)(Schatztruhe)[!][Apr 1999].iso / Aminet / gfx / misc / gnuplot-3.7src.lha / gnuplot-3.7src / gnuplot-3.7.lha / gnuplot-3.7 / hidden3d.c < prev    next >
Encoding:
C/C++ Source or Header  |  1998-12-09  |  84.8 KB  |  2,796 lines
  1. #ifndef lint
  2. static char *RCSid = "$Id: hidden3d.c,v 1.33 1998/06/18 14:55:10 ddenholm Exp $";
  3. #endif
  4.  
  5. /* GNUPLOT - hidden3d.c */
  6.  
  7. /*[
  8.  * Copyright 1986 - 1993, 1998   Thomas Williams, Colin Kelley
  9.  *
  10.  * Permission to use, copy, and distribute this software and its
  11.  * documentation for any purpose with or without fee is hereby granted,
  12.  * provided that the above copyright notice appear in all copies and
  13.  * that both that copyright notice and this permission notice appear
  14.  * in supporting documentation.
  15.  *
  16.  * Permission to modify the software is granted, but not the right to
  17.  * distribute the complete modified source code.  Modifications are to
  18.  * be distributed as patches to the released version.  Permission to
  19.  * distribute binaries produced by compiling modified sources is granted,
  20.  * provided you
  21.  *   1. distribute the corresponding source modifications from the
  22.  *    released version in the form of a patch file along with the binaries,
  23.  *   2. add special version identification to distinguish your version
  24.  *    in addition to the base release version number,
  25.  *   3. provide your name and address as the primary contact for the
  26.  *    support of your modified version, and
  27.  *   4. retain our contact information in regard to use of the base
  28.  *    software.
  29.  * Permission to distribute the released version of the source code along
  30.  * with corresponding source modifications in the form of a patch file is
  31.  * granted with same provisions 2 through 4 for binary distributions.
  32.  *
  33.  * This software is provided "as is" without express or implied warranty
  34.  * to the extent permitted by applicable law.
  35. ]*/
  36.  
  37.  
  38. /*
  39.  * 19 September 1992  Lawrence Crowl  (crowl@cs.orst.edu)
  40.  * Added user-specified bases for log scaling.
  41.  *
  42.  * 'someone'  contributed a complete rewrite of the hidden line
  43.  * stuff, for about beta 173 or so, called 'newhide.tgz'
  44.  *
  45.  * 1995-1997 Hans-Bernhard Br"oker (Broeker@physik.rwth-aachen.de)
  46.  *   Speedup, cleanup and several modification to integrate
  47.  *   'newhide' with 3.6 betas ~188 (start of work) up to 273.
  48.  *   As of beta 290, this is 'officially' in gnuplot.
  49.  *
  50.  */
  51.  
  52. #include "plot.h"
  53. #include "setshow.h"
  54.  
  55. /* TODO (HBB's notes, just in case you're interested):
  56.  * + fix all the problems I annotated by a 'FIXME' comment
  57.  * + Find out which value EPSILON should have, and why
  58.  * + Ask gnuplot-beta for a suitable way of concentrating
  59.  *   all those 'VERYSMALL', 'EPSILON', and other numbers
  60.  *   into one, consistent scheme, possibly based on
  61.  *   <float.h>. E.g., I'd say that for most applications,
  62.  *   the best 'epsilon' is either DBL_EPS or its square root.
  63.  * + redo all the profiling, esp. to find out if TEST_GRIDCHECK = 1
  64.  *   actually gains some speed, and if it is worth the memory
  65.  *   spent to store the bit masks
  66.  *   -> seems not improve speed at all, at least for my standard
  67.  *   test case (mandd.gpl) it even slows it down by ~ 10%
  68.  * + Evaluate the speed/memory comparison for storing vertex
  69.  *   indices instead of (double) floating point constants
  70.  *   to store {x|y}{min|max}
  71.  * + remove those of my comments that are now pointless
  72.  * + indent all that code
  73.  * + try to really understand that 'hl_buffer' stuff...
  74.  * + get rid of hardcoded things like sizeof(short) == 4,
  75.  *   those 0x7fff terms and all that.
  76.  * + restructure the hl_buffer storing method to make it more efficient
  77.  * + Try to cut the speed decrease of this code rel. to the old hidden
  78.  *   hidden line removal. (For 'mandd.gpl', it costs an overall
  79.  *   factor of 9 compared with my older versions of gnuplot!)
  80.  */
  81.  
  82. /*************************/
  83. /* Configuration section */
  84. /*************************/
  85.  
  86. /* HBB: if this module is compiled with TEST_GRIDCHECK = 1 defined,
  87.  * it will store the information about {x|y}{min|max} in an
  88.  * other (additional!) form: a bit mask, with each bit representing
  89.  * one horizontal or vertical strip of the screen. The bits
  90.  * for  strips a polygon spans are set to one. This allows to
  91.  * test for xy_overlap simply by comparing bit patterns.
  92.  */
  93. #ifndef TEST_GRIDCHECK
  94. #define TEST_GRIDCHECK 0
  95. #endif
  96.  
  97. /* HBB 961124; If you don't want the color-distinction between the
  98.  * 'top' and 'bottom' sides of the surface, like I do, then just compile
  99.  * with -DBACKSIDE_LINETYPE_OFFSET = 0. */
  100. #ifndef BACKSIDE_LINETYPE_OFFSET
  101. #define BACKSIDE_LINETYPE_OFFSET 1
  102. #endif
  103.  
  104. /* HBB 961127: defining FOUR_TRIANGLES = 1 will separate each quadrangle
  105.  * of data points into *four* triangles, by introducing an additional
  106.  * point at the mean of the four corners. Status: experimental 
  107.  */
  108. #ifndef FOUR_TRIANGLES
  109. #define FOUR_TRIANGLES 0
  110. #endif
  111.  
  112. /* HBB 970322: I removed the SENTINEL flag and its implementation
  113.  * completely.  Its speed improvement wasn't that great, and it never
  114.  * fully worked anyway, so that shouldn't make much of a difference at
  115.  * all. */
  116.  
  117. /* HBB 961212: this #define lets you choose if the diagonals that
  118.  * divide each original quadrangle in two triangles will be drawn
  119.  * visible or not: do draw them, define it to be 7L, otherwise let be
  120.  * 3L (for the FOUR_TRIANGLES method, the values are 7L and 1L) */
  121. /* drd : default to 3, for back-compatibility with 3.5 */
  122. #ifndef TRIANGLE_LINESDRAWN_PATTERN
  123. #define TRIANGLE_LINESDRAWN_PATTERN 3L
  124. #endif
  125.  
  126. /* HBB 970131: with HANDLE_UNDEFINED_POINTS = 1, let's try to handle
  127.  * UNDEFINED data points in the input we're given by the rest of
  128.  * gnuplot. We do that by marking these points by giving them z = -2
  129.  * (which normally never happens), and then refusing to build any
  130.  * polygon if one of its vertices has this mark. Thus, there's now a
  131.  * hole in the generated surface. */
  132. /* HBB 970322: some of this code is now hardwired (in
  133.  * PREPARE_POLYGON), so HANDLE_UNDEFINED_POINTS = 0 might have stopped
  134.  * to be a useful setting. I doubt anyone will miss it, anyway. */
  135. /* drd : 2 means undefined only. 1 means outrange treated as undefined */
  136. #ifndef HANDLE_UNDEFINED_POINTS
  137. #define HANDLE_UNDEFINED_POINTS 1
  138. #endif
  139. /* Symbolic value for 'do not handle Undefined Points specially' */
  140. #define UNHANDLED (UNDEFINED+1)
  141.  
  142. /* HBB 970322: If both subtriangles of a quad were cancelled, try if
  143.  * using the other diagonal is better. This only makes a difference if
  144.  * exactly one vertex of the quad is unusable, and that one is on the
  145.  * first tried diagonal. In such a case, the border of the hole in the
  146.  * surface will be less rough than with the previous method. To
  147.  * disable this, #define SHOW_ALTERNATIVE_DIAGONAL as 0 */
  148. #ifndef SHOW_ALTERNATIVE_DIAGONAL
  149. #define SHOW_ALTERNATIVE_DIAGONAL 1
  150. #endif
  151.  
  152. /* HBB 9790522: If the two triangles in a quad are both drawn, and
  153.  * they show different sides to the user (the quad is 'bent over'),
  154.  * then it often looks more sensible to draw the diagonal visibly to
  155.  * avoid other parts of the scene being obscured by what the user
  156.  * can't see. To disable, #define HANDLE_BENTOVER_QUADRANGLES 0 */
  157. #ifndef HANDLE_BENTOVER_QUADRANGLES
  158. #define HANDLE_BENTOVER_QUADRANGLES 1
  159. #endif
  160.  
  161. /* HBB 970618: The code of split_polygon_by_plane used to make a
  162.  * separate decision about what should happen to points that are 'in'
  163.  * the splitting plane (within EPSILON accuracy, i.e.), based on the
  164.  * orientation of the splitting plane. I had disabled that code for I
  165.  * assumed it might be the cause of some of the buggyness. I'm not yet
  166.  * fully convinced on wether that assumption holds or not, so it's now
  167.  * choosable.  OLD_SPLIT_PLANE == 1 will enable it. Some older comments
  168.  * from the source:*/
  169. /* HBB 970325: re-inserted this from older versions. Finally, someone
  170.  * came up with a test case where it is actually needed :-( */
  171. /* HBB 970427: seems to have been an incorrect analysis of that error.
  172.  * the problematic plot works without this as well. */
  173. #ifndef OLD_SPLIT_INPLANE
  174. #define OLD_SPLIT_INPLANE 1
  175. #endif
  176.  
  177. /* HBB 970618: this way, new edges introduced by splitting a polygon
  178.  * by the plane of another one will be made visible. Not too useful on
  179.  * production output, but can help in debugging. */
  180. #ifndef SHOW_SPLITTING_EDGES
  181. #define SHOW_SPLITTING_EDGES 0
  182. #endif
  183.  
  184. /********************************/
  185. /* END of configuration section */
  186. /********************************/
  187.  
  188.  
  189.  
  190. /* Variables to hold configurable values. This is meant to prepare for
  191.  * making these settable by the user via 'set hidden [option]...' */
  192.  
  193. int hiddenBacksideLinetypeOffset = BACKSIDE_LINETYPE_OFFSET;
  194. long hiddenTriangleLinesdrawnPattern = TRIANGLE_LINESDRAWN_PATTERN;
  195. int hiddenHandleUndefinedPoints = HANDLE_UNDEFINED_POINTS;
  196. int hiddenShowAlternativeDiagonal = SHOW_ALTERNATIVE_DIAGONAL;
  197. int hiddenHandleBentoverQuadrangles = HANDLE_BENTOVER_QUADRANGLES;
  198.  
  199. /* The functions to map from user 3D space into normalized -1..1 */
  200. #define map_x3d(x) ((x-min_array[FIRST_X_AXIS])*xscale3d-1.0)
  201. #define map_y3d(y) ((y-min_array[FIRST_Y_AXIS])*yscale3d-1.0)
  202. #define map_z3d(z) ((z-base_z)*zscale3d-1.0)
  203.  
  204. extern int suppressMove;
  205. extern int xright, xleft, ybot, ytop;
  206. extern int xmiddle, ymiddle, xscaler, yscaler;
  207. extern double min_array[], max_array[];
  208. extern int auto_array[], log_array[];
  209. extern double base_array[], log_base_array[];
  210. extern double xscale3d, yscale3d, zscale3d;
  211. extern double base_z;
  212.  
  213. extern int hidden_no_update, hidden_active;
  214. extern int hidden_line_type_above, hidden_line_type_below;
  215.  
  216. extern double trans_mat[4][4];
  217.  
  218.  
  219. #ifdef HAVE_STRINGIZE
  220. /* ANSI preprocessor concatenation */
  221. # define CONCAT(x,y) x##y
  222. # define CONCAT3(x,y,z) x##y##z
  223. #else
  224. /* K&R preprocessor concatenation */
  225. # define CONCAT(x,y) x/**/y
  226. # define CONCAT3(x,y,z) x/**/y/**/z
  227. #endif
  228.  
  229. /* Bitmap of the screen.  The array for each x value is malloc-ed as needed */
  230. /* HBB 961111: started parametrisation of type t_pnt, to prepare change from
  231.  * short to normal ints in **pnt. The other changes aren't always annotated,
  232.  * so watch out! */
  233. typedef unsigned short int t_pnt;
  234. #define PNT_BITS (CHAR_BIT*sizeof(t_pnt))
  235. /* caution! ANSI-dependant ! ? */
  236. #define PNT_MAX USHRT_MAX
  237. typedef t_pnt *tp_pnt;
  238. static tp_pnt *pnt;
  239.  
  240. /* Testroutine for the bitmap */
  241. /* HBB 961110: fixed slight problem indicated by lclint: 
  242.  * calling IS_UNSET with pnt == 0 (possible?) */
  243. /* HBB 961111: changed for new t_pnt, let compiler take care of optimisation
  244.  * `%' -> `&' */
  245. /* HBB 961124: switched semantics: as this was *always* called as !IS_SET(),
  246.  * it's better to check for *unset* bits right away: */
  247. #define IS_UNSET(X,Y) ((!pnt || pnt[X] == 0) ? 1 : !(((pnt[X])[(Y)/PNT_BITS] >> ((Y)%PNT_BITS)) & 0x01))
  248.  
  249. /* Amount of space we need for one vertical row of bitmap, 
  250.    and byte offset of first used element */
  251. static unsigned long y_malloc;
  252.  
  253. /* These numbers are chosen as dividers into the bitmap. */
  254. static int xfact, yfact;
  255. #define XREDUCE(X) ((X)/xfact)
  256. #define YREDUCE(Y) ((Y)/yfact)
  257.  
  258. /* These variables are only used for drawing the individual polygons
  259.    that make up the 3d figure.  After each polygon is drawn, the
  260.    information is copied to the bitmap: xmin_hl, ymin_hl are used to
  261.    keep track of the range of x values.  The arrays ymin_hl[],
  262.    ymax_hl[] are used to keep track of the minimum and maximum y
  263.    values used for each X value. */
  264.  
  265. /* HBB 961124: new macro, to avoid a wordsize-depence */
  266. #define HL_EXTENT_X_MAX UINT_MAX    /* caution! ANSI-only !? */
  267. static unsigned int xmin_hl, xmax_hl;
  268.  
  269. /* HBB: parametrized type of hl_buffer elements, to allow easier
  270.    changing later on: */
  271. #define HL_EXTENT_Y_MAX SHRT_MAX    /* caution! ANSI-only !? */
  272. typedef short int t_hl_extent_y;
  273. static t_hl_extent_y *ymin_hl, *ymax_hl;
  274.  
  275.  
  276. /* hl_buffer is a buffer which is needed to draw polygons with very small 
  277.    angles correct:
  278.  
  279.    One polygon may be split during the sorting algorithmus into
  280.    several polygons. Before drawing a part of a polygon, I save in
  281.    this buffer all lines of the polygon, which will not yet be drawn.
  282.    If then the actual part draws a virtual line (a invisible line,
  283.    which splits the polygon into 2 parts), it draws the line visible
  284.    in those points, which are set in the buffer.
  285.  
  286.    The layout of the buffer:
  287.    hl_buffer[i] stands for the i'th column of the bitmap.
  288.    Each column is splitted into pairs of points (a,b).
  289.    Each pair describes a cross of one line with the column. */
  290.  
  291. struct Cross {
  292.     int a, b;
  293.     struct Cross GPHUGE *next;
  294. };
  295.  
  296. static struct Cross GPHUGE *GPHUGE * hl_buffer;
  297.  
  298. /* HBB 980303: added a global array of Cross structures, to save lots
  299.  * of gp_alloc() calls (3 millions in a run through all.dem!)  */
  300.  
  301. #define CROSS_STORE_INCREASE 500    /* number of Cross'es to alloc at a time */
  302. static struct Cross *Cross_store = 0;
  303. static int last_Cross_store = 0, free_Cross_store = 0;
  304.  
  305. struct Vertex {
  306.     coordval x, y, z;
  307.     int style;
  308. };
  309.  
  310. /* HBB 971114: two new macros, to properly hide the machinery that
  311.  *implements flagging of vertices as 'undefined' (or 'out of range',
  312.  *handled equally). Thanks to Lars Hecking for pointing this out */
  313. #define FLAG_VERTEX_AS_UNDEFINED(v) \
  314.   do { (v).z = -2.0; } while (0)
  315. #define VERTEX_IS_UNDEFINED(v) ((v).z == -2.0)
  316.  
  317. /* These values serve as flags to detect loops of polygons obscuring
  318.  * each other in in_front() */
  319. typedef enum {
  320.     is_untested = 0, is_tested, is_in_loop, is_moved_or_split
  321. } t_poly_tested;
  322.  
  323. /* Strings for printing out values of type t_poly_tested: */
  324. static const char *string_poly_tested[] =
  325. {
  326.     "is_untested",
  327.     "is_tested",
  328.     "is_in_loop",
  329.     "is_moved_or_split"
  330. };
  331.  
  332. struct Polygon {
  333.     int n;            /* amount of vertices */
  334.     long *vertex;        /* The vertices (indices on vlist) */
  335.     long line;            /* i'th bit shows, if the i'th line should be drawn */
  336.     coordval xmin, xmax, ymin, ymax, zmin, zmax;
  337.     /* min/max in all three directions */
  338.     struct lp_style_type *lp;    /* pointer to l/p properties */
  339.     int style;            /* The style of the lines */
  340.     long next;            /* index of next polygon in z-sorted list */
  341.     long next_frag;        /* next fragment of same polygon... */
  342.     int id;            /* Which polygons belong together? */
  343.     t_poly_tested tested;    /* To determine a loop during the sorting algorithm */
  344.     /* HBB 980317: on the 16 bit PC platforms, the struct size must
  345.      * be a power of two, so it exactly fits into a 64KB segment. First
  346.      * we'll add the TEST_GRIDCHECK fields, regardless wether that
  347.      * feature was activated or not. */
  348. #if TEST_GRIDCHECK || defined(DOS16) || defined(WIN16)
  349.     unsigned int xextent, yextent;
  350. #endif
  351.     /* HBB 980317: the remaining bit of padding. */
  352. #if defined(DOS16) || defined(WIN16) || defined(WIN32)
  353.     char dummies[8];
  354. #endif
  355. };
  356.  
  357. typedef enum {            /* result type for polygon_plane_intersection() */
  358.     infront, inside, behind, intersects
  359. } t_poly_plane_intersect;
  360.  
  361. static struct Vertex GPHUGE *vlist;    /* The vertices */
  362. static struct Polygon GPHUGE *plist;    /* The polygons */
  363. static long nvert, npoly;    /* amount of vertices and polygons */
  364. static long pfree, vert_free;    /* index on the first free entries */
  365.  
  366. static long pfirst;        /* Index on the first polygon */
  367.  
  368. /* macros for (somewhat) safer handling of the dynamic array of
  369.  * polygons, 'plist[]' */
  370. #define EXTEND_PLIST() \
  371.     plist = (struct Polygon GPHUGE *) gp_realloc(plist, \
  372.       (unsigned long)sizeof(struct Polygon)*(npoly += 10L), "hidden plist")
  373.  
  374. #define CHECK_PLIST() if (pfree >= npoly) EXTEND_PLIST()
  375. #define CHECK_PLIST_2() if (pfree+1 >= npoly) EXTEND_PLIST()
  376. #define CHECK_PLIST_EXTRA(extra) if (pfree >= npoly) { EXTEND_PLIST() ; extra; }
  377.  
  378.  
  379. /* precision of calculations in normalized space: */
  380. #define EPSILON 1e-5        /* HBB: was 1.0E-5 */
  381. #define SIGNOF(X)  ( ((X)<-EPSILON) ? -1: ((X)>EPSILON) )
  382.  
  383. /* Some inexact relations: == , > , >= */
  384. #define EQ(X,Y)  (fabs( (X) - (Y) ) < EPSILON)    /* X == Y */
  385. #define GR(X,Y)  ((X) >  (Y) + EPSILON)        /* X >  Y */
  386. #define GE(X,Y)  ((X) >= (Y) - EPSILON)        /* X >= Y */
  387.  
  388. /* Test, if two vertices are near each other */
  389. #define V_EQUAL(a,b) ( GE(0.0, fabs((a)->x - (b)->x) + \
  390.    fabs((a)->y - (b)->y) + fabs((a)->z - (b)->z)) )
  391.  
  392. #define SETBIT(a,i,set) a= (set) ?  (a) | (1L<<(i)) : (a) & ~(1L<<(i))
  393. #define GETBIT(a,i) (((a)>>(i))&1L)
  394.  
  395. #define UINT_BITS (CHAR_BIT*sizeof(unsigned int))
  396. #define USHRT_BITS (CHAR_BIT*sizeof(unsigned short))
  397.  
  398.  
  399. static void print_polygon __PROTO((struct Polygon GPHUGE * p, const char *pname));
  400.  
  401. /* HBB: new routine, to help debug 'Logic errors', mainly */
  402. static void print_polygon(p, pname)
  403. struct Polygon GPHUGE *p;
  404. const char *pname;
  405. {
  406.     struct Vertex GPHUGE *v;
  407.     int i;
  408.  
  409.     fprintf(stderr, "#%s:(ind %d) n:%d, id:%d, next:%ld, tested:%s\n",
  410.         pname, (int) (p - plist), p->n, p->id, p->next,
  411.         string_poly_tested[(int) (p->tested)]);
  412.     fprintf(stderr, "#zmin:%g, zmax:%g\n", p->zmin, p->zmax);
  413.     for (i = 0; i < p->n; i++, v++) {
  414.     v = vlist + p->vertex[i];
  415.     fprintf(stderr, "%g %g %g \t%ld\n", v->x, v->y, v->z, p->vertex[i]);
  416.     }
  417.     /*repeat first vertex, so the line gets closed: */
  418.     v = vlist + p->vertex[0];
  419.     fprintf(stderr, "%g %g %g \t%ld\n", v->x, v->y, v->z, p->vertex[0]);
  420.     /*two blank lines, as a multimesh separator: */
  421.     fputs("\n\n", stderr);
  422. }
  423.  
  424. /* Gets Minimum 'C' value of polygon, C is x, y, or z: */
  425. #define GET_MIN(p, C, min) do { \
  426.   int i; \
  427.   long *v = p->vertex; \
  428.                        \
  429.   min = vlist[*v++].C; \
  430.   for (i = 1; i< p->n; i++, v++) \
  431.     if (vlist[*v].C < min) \
  432.       min = vlist[*v].C; \
  433. } while (0)
  434.  
  435. /* Gets Maximum 'C' value of polygon, C is x, y, or z: */
  436. #define GET_MAX(p, C, max) do { \
  437.   int i; \
  438.   long *v = p->vertex; \
  439.                       \
  440.   max = vlist[*v++].C; \
  441.   for (i = 1; i< p->n; i++, v++) \
  442.     if (vlist[*v].C > max) \
  443.       max = vlist[*v].C; \
  444. } while (0)
  445.  
  446. /* Prototypes for internal functions of this module. */
  447. static void map3d_xyz __PROTO((double x, double y, double z,
  448.                    struct Vertex GPHUGE * v));
  449. static int reduce_polygon __PROTO((int *n, long **v, long *line, int nmin));
  450. static void build_polygon __PROTO((struct Polygon GPHUGE * p, int n,
  451.         long *v, long line, int style, struct lp_style_type * lp,
  452.            long next, long next_frag, int id, t_poly_tested tested));
  453. static GP_INLINE int maybe_build_polygon __PROTO((struct Polygon GPHUGE * p,
  454.      int n, long *v, long line, int style, struct lp_style_type * lp,
  455.            long next, long next_frag, int id, t_poly_tested tested));
  456. static void init_polygons __PROTO((struct surface_points * plots, int pcount));
  457. static int compare_by_zmax __PROTO((const void *p1, const void *p2));
  458. static void sort_by_zmax __PROTO((void));
  459. static int obscured __PROTO((struct Polygon GPHUGE * p));
  460. static GP_INLINE int xy_overlap __PROTO((struct Polygon GPHUGE * a, struct Polygon GPHUGE * b));
  461. static void get_plane __PROTO((struct Polygon GPHUGE * p, double *a, double *b,
  462.                    double *c, double *d));
  463. static t_poly_plane_intersect polygon_plane_intersection __PROTO((struct Polygon GPHUGE * p, double a, double b, double c, double d));
  464. static int intersect __PROTO((struct Vertex GPHUGE * v1,
  465.                   struct Vertex GPHUGE * v2, struct Vertex GPHUGE * v3, struct Vertex GPHUGE * v4));
  466. static int v_intersect __PROTO((struct Vertex GPHUGE * v1,
  467.                 struct Vertex GPHUGE * v2, struct Vertex GPHUGE * v3, struct Vertex GPHUGE * v4));
  468. static int intersect_polygon __PROTO((struct Vertex GPHUGE * v1, struct Vertex GPHUGE * v2, struct Polygon GPHUGE * p));
  469. static int full_xy_overlap __PROTO((struct Polygon GPHUGE * a, struct Polygon GPHUGE * b));
  470. static long build_new_vertex __PROTO((long V1, long V2, double w));
  471. static long split_polygon_by_plane __PROTO((long P, double a, double b,
  472.                     double c, double d, TBOOLEAN f));
  473. static int in_front __PROTO((long Last, long Test));
  474.  
  475. /* HBB 980303: new, for the new back-buffer for *Cross structures: */
  476. static GP_INLINE struct Cross *get_Cross_from_store __PROTO((void));
  477. static GP_INLINE void init_Cross_store __PROTO((void));
  478.  
  479. static GP_INLINE TBOOLEAN hl_buffer_set __PROTO((int xv, int yv));
  480. static GP_INLINE void update_hl_buffer_column __PROTO((int xv, int ya, int yv));
  481. static void draw_clip_line_buffer __PROTO((int x1, int y1, int x2, int y2));
  482. static void draw_clip_line_update __PROTO((int x1, int y1, int x2, int y2,
  483.                        int virtual));
  484. static GP_INLINE void clip_vector_h __PROTO((int x, int y));
  485. static GP_INLINE void clip_vector_virtual __PROTO((int x, int y));
  486. static GP_INLINE void clip_vector_buffer __PROTO((int x, int y));
  487. static void draw __PROTO((struct Polygon GPHUGE * p));
  488.  
  489.  
  490. /* Set the options for hidden3d. To be called from set.c, when the
  491.  * user has begun a command with 'set hidden3d', to parse the rest of
  492.  * that command */
  493. void set_hidden3doptions()
  494. {
  495.     struct value t;
  496.  
  497.     c_token++;
  498.  
  499.     if (END_OF_COMMAND) {
  500.     return;
  501.     }
  502.     if (almost_equals(c_token, "def$aults")) {
  503.     /* reset all parameters to defaults */
  504.     hiddenBacksideLinetypeOffset = BACKSIDE_LINETYPE_OFFSET;
  505.     hiddenTriangleLinesdrawnPattern = TRIANGLE_LINESDRAWN_PATTERN;
  506.     hiddenHandleUndefinedPoints = HANDLE_UNDEFINED_POINTS;
  507.     hiddenShowAlternativeDiagonal = SHOW_ALTERNATIVE_DIAGONAL;
  508.     hiddenHandleBentoverQuadrangles = HANDLE_BENTOVER_QUADRANGLES;
  509.     c_token++;
  510.     if (!END_OF_COMMAND)
  511.         int_error("No further options allowed after 'defaults'", c_token);
  512.     return;
  513.     }
  514.     if (almost_equals(c_token, "off$set")) {
  515.     c_token++;
  516.     hiddenBacksideLinetypeOffset = (int) real(const_express(&t));
  517.     } else if (almost_equals(c_token, "nooff$set")) {
  518.     c_token++;
  519.     hiddenBacksideLinetypeOffset = 0;
  520.     }
  521.     if (END_OF_COMMAND) {
  522.     return;
  523.     }
  524.     if (almost_equals(c_token, "tri$anglepattern")) {
  525.     c_token++;
  526.     hiddenTriangleLinesdrawnPattern = (int) real(const_express(&t));
  527.     }
  528.     if (END_OF_COMMAND) {
  529.     return;
  530.     }
  531.     if (almost_equals(c_token, "undef$ined")) {
  532.     int tmp;
  533.  
  534.     c_token++;
  535.     tmp = (int) real(const_express(&t));
  536.     if (tmp <= 0 || tmp > UNHANDLED)
  537.         tmp = UNHANDLED;
  538.     hiddenHandleUndefinedPoints = tmp;
  539.     } else if (almost_equals(c_token, "nound$efined")) {
  540.     c_token++;
  541.     hiddenHandleUndefinedPoints = UNHANDLED;
  542.     }
  543.     if (END_OF_COMMAND) {
  544.     return;
  545.     }
  546.     if (almost_equals(c_token, "alt$diagonal")) {
  547.     c_token++;
  548.     hiddenShowAlternativeDiagonal = 1;
  549.     } else if (almost_equals(c_token, "noalt$diagonal")) {
  550.     c_token++;
  551.     hiddenShowAlternativeDiagonal = 0;
  552.     }
  553.     if (END_OF_COMMAND) {
  554.     return;
  555.     }
  556.     if (almost_equals(c_token, "bent$over")) {
  557.     c_token++;
  558.     hiddenHandleBentoverQuadrangles = 1;
  559.     } else if (almost_equals(c_token, "nobent$over")) {
  560.     c_token++;
  561.     hiddenHandleBentoverQuadrangles = 0;
  562.     }
  563.     if (!END_OF_COMMAND) {
  564.     int_error("No such option to hidden3d (or wrong order)", c_token);
  565.     }
  566. }
  567.  
  568. /* HBB 970619: new routine for displaying the option settings */
  569. void show_hidden3doptions()
  570. {
  571.     fprintf(stderr,"\
  572. \t  Back side of surfaces has linestyle offset of %d\n\
  573. \t  Bit-Mask of Lines to draw in each triangle is %ld\n\
  574. \t  %d: ",
  575.         hiddenBacksideLinetypeOffset, hiddenTriangleLinesdrawnPattern,
  576.         hiddenHandleUndefinedPoints);
  577.  
  578.     switch (hiddenHandleUndefinedPoints) {
  579.     case OUTRANGE:
  580.     fputs("Outranged and undefined datapoints are omitted from the surface.\n", stderr);
  581.     break;
  582.     case UNDEFINED:
  583.     fputs("Only undefined datapoints are omitted from the surface.\n", stderr);
  584.     break;
  585.     case UNHANDLED:
  586.     fputs("Will not check for undefined datapoints (may cause crashes).\n", stderr);
  587.     break;
  588.     default:
  589.     fputs("This value is illegal!!!\n", stderr);
  590.     break;
  591.     }
  592.     fprintf(stderr,"\
  593. \t  Will %suse other diagonal if it gives a less jaggy outline\n\
  594. \t  Will %sdraw diagonal visibly if quadrangle is 'bent over'\n",
  595.         hiddenShowAlternativeDiagonal ? "" : "not ",
  596.         hiddenHandleBentoverQuadrangles ? "" : "not ");
  597. }
  598.  
  599. /* HBB 971117: new function. */
  600. /* Implements proper 'save'ing of the new hidden3d options... */
  601. void save_hidden3doptions(fp)
  602. FILE *fp;
  603. {
  604.     if (!hidden3d) {
  605.     fputs("set nohidden3d\n", fp);
  606.     return;
  607.     }
  608.     fprintf(fp, "set hidden3d offset %d trianglepattern %ld undefined %d %saltdiagonal %sbentover\n",
  609.         hiddenBacksideLinetypeOffset,
  610.         hiddenTriangleLinesdrawnPattern,
  611.         hiddenHandleUndefinedPoints,
  612.         hiddenShowAlternativeDiagonal ? "" : "no",
  613.         hiddenHandleBentoverQuadrangles ? "" : "no");
  614. }
  615.  
  616. /* Initialize the necessary steps for hidden line removal and
  617.    initialize global variables. */
  618. void init_hidden_line_removal()
  619. {
  620.     int i;
  621.  
  622.     /* Check for some necessary conditions to be set elsewhere: */
  623.     /* HandleUndefinedPoints mechanism depends on these: */
  624.     assert(OUTRANGE == 1);
  625.     assert(UNDEFINED == 2);
  626.     /* '3' makes the test easier in the critical section */
  627.     if (hiddenHandleUndefinedPoints < OUTRANGE)
  628.     hiddenHandleUndefinedPoints = UNHANDLED;
  629.  
  630.     /*  We want to keep the bitmap size less than 2048x2048, so we choose
  631.      *  integer dividers for the x and y coordinates to keep the x and y
  632.      *  ranges less than 2048.  In practice, the x and y sizes for the bitmap
  633.      *  will be somewhere between 1024 and 2048, except in cases where the
  634.      *  coordinates ranges for the device are already less than 1024.
  635.      *  We do this mainly to control the size of the bitmap, but it also
  636.      *  speeds up the computation.  We maintain separate dividers for
  637.      *  x and y.
  638.      */
  639.     xfact = (xright - xleft) / 1024;
  640.     yfact = (ytop - ybot) / 1024;
  641.     if (xfact == 0)
  642.     xfact = 1;
  643.     if (yfact == 0)
  644.     yfact = 1;
  645.     if (pnt == 0) {
  646.     i = XREDUCE(xright) - XREDUCE(xleft) + 1;
  647.     pnt = (tp_pnt *) gp_alloc(
  648.              (unsigned long) (i * sizeof(tp_pnt)), "hidden pnt");
  649.     while (--i >= 0)
  650.         pnt[i] = (tp_pnt) 0;
  651.     }
  652. }
  653.  
  654. /* Reset the hidden line data to a fresh start.                              */
  655. void reset_hidden_line_removal()
  656. {
  657.     int i;
  658.     if (pnt) {
  659.     for (i = 0; i <= XREDUCE(xright) - XREDUCE(xleft); i++) {
  660.         if (pnt[i]) {
  661.         free(pnt[i]);
  662.         pnt[i] = 0;
  663.         };
  664.     };
  665.     };
  666. }
  667.  
  668.  
  669. /* Terminates the hidden line removal process.                  */
  670. /* Free any memory allocated by init_hidden_line_removal above. */
  671. void term_hidden_line_removal()
  672. {
  673.     if (pnt) {
  674.     int j;
  675.     for (j = 0; j <= XREDUCE(xright) - XREDUCE(xleft); j++) {
  676.         if (pnt[j]) {
  677.         free(pnt[j]);
  678.         pnt[j] = 0;
  679.         };
  680.     };
  681.     free(pnt);
  682.     pnt = 0;
  683.     };
  684.     if (ymin_hl)
  685.     free(ymin_hl), ymin_hl = 0;
  686.     if (ymax_hl)
  687.     free(ymax_hl), ymax_hl = 0;
  688. }
  689.  
  690.  
  691. /* Maps from normalized space to terminal coordinates */
  692. #define TERMCOORD(v,x,y) x = ((int)((v)->x * xscaler)) + xmiddle, \
  693.              y = ((int)((v)->y * yscaler)) + ymiddle;
  694.  
  695.  
  696. static void map3d_xyz(x, y, z, v)
  697. double x, y, z;            /* user coordinates */
  698. struct Vertex GPHUGE *v;    /* the point in normalized space */
  699. {
  700.     int i, j;
  701.     double V[4], Res[4],    /* Homogeneous coords. vectors. */
  702.      w = trans_mat[3][3];
  703.     /* Normalize object space to -1..1 */
  704.     V[0] = map_x3d(x);
  705.     V[1] = map_y3d(y);
  706.     V[2] = map_z3d(z);
  707.     V[3] = 1.0;
  708.     Res[0] = 0;            /*HBB 961110: lclint wanted this... Why? */
  709.     for (i = 0; i < 3; i++) {
  710.     Res[i] = trans_mat[3][i];    /* Initiate it with the weight factor. */
  711.     for (j = 0; j < 3; j++)
  712.         Res[i] += V[j] * trans_mat[j][i];
  713.     }
  714.     for (i = 0; i < 3; i++)
  715.     w += V[i] * trans_mat[i][3];
  716.     if (w == 0)
  717.     w = 1.0e-5;
  718.     v->x = Res[0] / w;
  719.     v->y = Res[1] / w;
  720.     v->z = Res[2] / w;
  721. }
  722.  
  723.  
  724. /* remove unnecessary Vertices from a polygon: */
  725. static int reduce_polygon(n, v, line, nmin)
  726. int *n;                /* number of vertices */
  727. long **v;            /* the vertices */
  728. long *line;            /* information which line should be drawn */
  729. int nmin;            /* if the number reduces under nmin, forget polygon */
  730.      /* Return value 1: the reduced polygon has still more or equal 
  731.         vertices than nmin.
  732.         Return value 0: the polygon was trivial; memory is given free */
  733. {
  734.     int less, i;
  735.     long *w = *v;
  736.     for (less = 0, i = 0; i < *n - 1; i++) {
  737.     if (V_EQUAL(vlist + w[i], vlist + w[i + 1])
  738.     /* HBB 970321: try to remove an endless loop bug (part1/2) */
  739.         && ((w[i] == w[i + 1])
  740.         || !GETBIT(*line, i)
  741.         || GETBIT(*line, i + 1)
  742.         || ((i > 0) ? GETBIT(*line, i - 1) : GETBIT(*line, *n - 1))))
  743.         less++;
  744.     else if (less) {
  745.         w[i - less] = w[i];
  746.         SETBIT(*line, i - less, GETBIT(*line, i));
  747.     }
  748.     }
  749.     if (i - less > 0
  750.     && V_EQUAL(vlist + w[i], vlist + w[0])
  751.     /* HBB 970321: try to remove an endless loop bug (part2/2) */
  752.     && (w[i] == w[0]
  753.         || !GETBIT(*line, i)
  754.         || GETBIT(*line, i - 1)
  755.         || GETBIT(*line, 0)))
  756.     less++;
  757.     else if (less) {
  758.     w[i - less] = w[i];
  759.     SETBIT(*line, i - less, GETBIT(*line, i));
  760.     }
  761.     *n -= less;
  762.     if (*n < nmin) {
  763.     free(w);
  764.     *v = 0;            /* HBB 980213: signal that v(==w) is free()d */
  765.     return 0;
  766.     }
  767.     if (less) {
  768.     w = (long *) gp_realloc(w, (unsigned long) sizeof(long) * (*n),
  769.                 "hidden red_poly");
  770.     *v = w;
  771.     for (; less > 0; less--)
  772.         SETBIT(*line, *n + less - 1, 0);
  773.     }
  774.     return 1;
  775. }
  776.  
  777. /* Write polygon in plist */
  778. /*
  779.  * HBB: changed this to precalculate {x|y}{min|max}
  780.  */
  781. static void build_polygon(p, n, v, line, style, lp, next, next_frag, id, tested)
  782. struct Polygon GPHUGE *p;    /* this should point at a free entry in plist */
  783. int n;                /* number of vertices */
  784. long *v;            /* array of indices on the vertices (in vlist) */
  785. long line;            /* information, which line should be drawn */
  786. int style;            /* color of the lines */
  787. struct lp_style_type *lp;    /* pointer to line/point properties */
  788. long next;            /* next polygon in the list */
  789. long next_frag;            /* next fragment of same polygon */
  790. int id;                /* Which polygons belong together? */
  791. t_poly_tested tested;        /* Is polygon already on the right place of list? */
  792. {
  793.     int i;
  794.     struct Vertex vert;
  795.  
  796.     if (p >= plist + npoly)
  797.     fputs("uh oh !\n", stderr);
  798.  
  799.     CHECK_POINTER(plist, p);
  800.  
  801.     p->n = n;
  802.     p->vertex = v;
  803.     p->line = line;
  804.     p->lp = lp;            /* line and point properties */
  805.     p->style = style;
  806.     CHECK_POINTER(vlist, (vlist + v[0]));
  807.     vert = vlist[v[0]];
  808.     p->xmin = p->xmax = vert.x;
  809.     p->ymin = p->ymax = vert.y;
  810.     p->zmin = p->zmax = vert.z;
  811.     for (i = 1; i < n; i++) {
  812.     CHECK_POINTER(vlist, (vlist + v[i]));
  813.     vert = vlist[v[i]];
  814.     if (vert.x < p->xmin)
  815.         p->xmin = vert.x;
  816.     else if (vert.x > p->xmax)
  817.         p->xmax = vert.x;
  818.     if (vert.y < p->ymin)
  819.         p->ymin = vert.y;
  820.     else if (vert.y > p->ymax)
  821.         p->ymax = vert.y;
  822.     if (vert.z < p->zmin)
  823.         p->zmin = vert.z;
  824.     else if (vert.z > p->zmax)
  825.         p->zmax = vert.z;
  826.     }
  827.  
  828. #if TEST_GRIDCHECK
  829. /* FIXME: this should probably use a table of bit-patterns instead... */
  830.     p->xextent = (~(~0U << (unsigned int) ((p->xmax + 1.0) / 2.0 * UINT_BITS + 1.0)))
  831.     & (~0U << (unsigned int) ((p->xmin + 1.0) / 2.0 * UINT_BITS));
  832.     p->yextent = (~(~0U << (unsigned int) ((p->ymax + 1.0) / 2.0 * UINT_BITS + 1.0)))
  833.     & (~0U << (unsigned int) ((p->ymin + 1.0) / 2.0 * UINT_BITS));
  834. #endif
  835.  
  836.     p->next = next;
  837.     p->next_frag = next_frag;    /* HBB 961201: link fragments, to speed up draw() q-loop */
  838.     p->id = id;
  839.     p->tested = tested;
  840. }
  841.  
  842.  
  843. /* get the amount of curves in a plot */
  844. #define GETNCRV(NCRV) do {\
  845.    if(this_plot->plot_type == FUNC3D)        \
  846.      for(icrvs = this_plot->iso_crvs,NCRV = 0; \
  847.      icrvs;icrvs = icrvs->next,NCRV++);  \
  848.    else if(this_plot->plot_type == DATA3D) \
  849.         NCRV = this_plot->num_iso_read;    \
  850.     else {                                 \
  851.         graph_error("Plot type is neither function nor data"); \
  852.         return; /* stops a gcc -Wunitialised warning */        \
  853.     } \
  854. } while (0)
  855.  
  856. /* Do we see the top or bottom of the polygon, or is it 'on edge'? */
  857. #define GET_SIDE(vlst,csign) do { \
  858.   double c = vlist[vlst[0]].x * (vlist[vlst[1]].y - vlist[vlst[2]].y) + \
  859.     vlist[vlst[1]].x * (vlist[vlst[2]].y - vlist[vlst[0]].y) + \
  860.     vlist[vlst[2]].x * (vlist[vlst[0]].y - vlist[vlst[1]].y); \
  861.   csign = SIGNOF (c); \
  862. } while (0)
  863.  
  864. /* HBB 970131: new function, to allow handling of undefined datapoints
  865.  * by leaving out the corresponding polygons from the list.
  866.  * Return Value: 1 if polygon was built, 0 otherwise */
  867. static GP_INLINE int maybe_build_polygon(p, n, v, line, style, lp,
  868.                      next, next_frag, id, tested)
  869. struct Polygon GPHUGE *p;
  870. struct lp_style_type *lp;
  871. int n, style, id;
  872. t_poly_tested tested;
  873. long *v, line, next, next_frag;
  874. {
  875.     int i;
  876.  
  877.     assert(v);
  878.  
  879.     if (hiddenHandleUndefinedPoints != UNHANDLED)
  880.     for (i = 0; i < n; i++)
  881.         if (VERTEX_IS_UNDEFINED(vlist[v[i]])) {
  882.         free(v);    /* HBB 980213: free mem... */
  883.         v = 0;
  884.         return 0;    /* *don't* build the polygon! */
  885.         }
  886.     build_polygon(p, n, v, line, style, lp, next, next_frag, id, tested);
  887.     return 1;
  888. }
  889.  
  890. /* HBB 970322: new macro, invented because code like this occured several
  891.  * times inside init_polygons, esp. after I added the option to choose the
  892.  * other diagonal to divide a quad when necessary */
  893. /* HBB 980213: Here and below in init_polygons, added code to properly
  894.  * free any allocated vertex lists ('v' here), or avoid allocating
  895.  * them in the first place. Thanks to Stefan Schroepfer for reporting
  896.  * the leak.  */
  897. #define PREPARE_POLYGON(n,v,i0,i1,i2,line,c,border,i_chk,ncrv_chk,style) do {\
  898.   (n) = 3; \
  899.   if (VERTEX_IS_UNDEFINED(vlist[vert_free + (i0)])\
  900.       || VERTEX_IS_UNDEFINED(vlist[vert_free + (i1)]) \
  901.       || VERTEX_IS_UNDEFINED(vlist[vert_free + (i2)])) { \
  902.     /* These values avoid any further action on this polygon */\
  903.     (v) = 0; /* flag this as undefined */ \
  904.     (c) = (border) = 0; \
  905.   } else {\
  906.   (v) = (long *) gp_alloc ((unsigned long) sizeof (long) * (n), "hidden PREPARE_POLYGON"); \
  907.   (v)[0] = vert_free + (i0);\
  908.   (v)[1] = vert_free + (i1);\
  909.   (v)[2] = vert_free + (i2);\
  910.   (line) = hiddenTriangleLinesdrawnPattern;\
  911.     GET_SIDE((v),(c));\
  912.     /* Is this polygon at the border of the surface? */\
  913.     (border) = (i == (i_chk) || ncrv == (ncrv_chk));\
  914.   }\
  915.   (style) = ((c) >= 0) ? above : below;\
  916. } while (0)
  917.  
  918. static void init_polygons(plots, pcount)
  919. struct surface_points *plots;
  920. int pcount;
  921.      /* Initialize vlist, plist */
  922. {
  923.     long int i;
  924.     struct surface_points *this_plot;
  925.     int surface;
  926.     long int ncrv, ncrv1;
  927.     struct iso_curve *icrvs;
  928.     int row_offset;
  929.     int id;
  930.     int n1, n2;            /* number of vertices of a Polygon */
  931.     long *v1, *v2;        /* the Vertices */
  932.     long line1, line2;        /* which  lines should be drawn */
  933.     int above = -1, below, style1, style2;    /* the line color */
  934.     struct lp_style_type *lp;    /* pointer to line and point properties */
  935.     int c1, c2;            /* do we see the front or the back */
  936.     TBOOLEAN border1, border2;    /* is it at the border of the surface */
  937.  
  938.     /* allocate memory for polylist and nodelist */
  939.     nvert = npoly = 0;
  940.     for (this_plot = plots, surface = 0;
  941.      surface < pcount;
  942.      this_plot = this_plot->next_sp, surface++) {
  943.     GETNCRV(ncrv);
  944.     switch (this_plot->plot_style) {
  945.     case LINESPOINTS:
  946.     case STEPS:        /* handle as LINES */
  947.     case FSTEPS:
  948.     case HISTEPS:
  949.     case LINES:
  950. #if FOUR_TRIANGLES
  951.         nvert += ncrv * (2 * this_plot->iso_crvs->p_count - 1);
  952. #else
  953.         nvert += ncrv * (this_plot->iso_crvs->p_count);
  954. #endif
  955.         npoly += 2 * (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  956.         break;
  957.     case DOTS:
  958.     case XERRORBARS:    /* handle as POINTSTYLE */
  959.     case YERRORBARS:
  960.     case XYERRORBARS:
  961.     case BOXXYERROR:
  962.     case BOXERROR:
  963.     case CANDLESTICKS:    /* HBB: these as well */
  964.     case FINANCEBARS:
  965.     case VECTOR:
  966.     case POINTSTYLE:
  967.         nvert += ncrv * (this_plot->iso_crvs->p_count);
  968.         npoly += (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  969.         break;
  970.     case BOXES:        /* handle as IMPULSES */
  971.     case IMPULSES:
  972.         nvert += 2 * ncrv * (this_plot->iso_crvs->p_count);
  973.         npoly += (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  974.         break;
  975.     }
  976.     }
  977.     vlist = (struct Vertex GPHUGE *)
  978.     gp_alloc((unsigned long) sizeof(struct Vertex) * nvert, "hidden vlist");
  979.     plist = (struct Polygon GPHUGE *)
  980.     gp_alloc((unsigned long) sizeof(struct Polygon) * npoly, "hidden plist");
  981.  
  982.     /* initialize vlist: */
  983.     for (vert_free = 0, this_plot = plots, surface = 0;
  984.      surface < pcount;
  985.      this_plot = this_plot->next_sp, surface++) {
  986.     switch (this_plot->plot_style) {
  987.     case LINESPOINTS:
  988.     case BOXERROR:        /* handle as POINTSTYLE */
  989.     case BOXXYERROR:
  990.     case XERRORBARS:
  991.     case YERRORBARS:
  992.     case XYERRORBARS:
  993.     case CANDLESTICKS:    /* HBB: these as well */
  994.     case FINANCEBARS:
  995.     case VECTOR:
  996.     case POINTSTYLE:
  997.         above = this_plot->lp_properties.p_type;
  998.         break;
  999.     case STEPS:        /* handle as LINES */
  1000.     case FSTEPS:
  1001.     case HISTEPS:
  1002.     case LINES:
  1003.     case DOTS:
  1004.     case BOXES:        /* handle as IMPULSES */
  1005.     case IMPULSES:
  1006.         above = -1;
  1007.         break;
  1008.     }
  1009.     GETNCRV(ncrv1);
  1010.     for (ncrv = 0, icrvs = this_plot->iso_crvs;
  1011.          ncrv < ncrv1 && icrvs;
  1012.          ncrv++, icrvs = icrvs->next) {
  1013.         struct coordinate GPHUGE *points = icrvs->points;
  1014.  
  1015.         for (i = 0; i < icrvs->p_count; i++) {
  1016.         /* As long as the point types keep OUTRANGE == 1 and
  1017.          * UNDEFINED == 2, we can just compare
  1018.          * hiddenHandleUndefinedPoints to the type of the point. To
  1019.          * simplify this, let UNHANDLED := UNDEFINED+1 mean 'do not
  1020.          * handle undefined or outranged points at all' (dangerous
  1021.          * option, that is...) */
  1022.         if (points[i].type >= hiddenHandleUndefinedPoints) {
  1023.             /* mark this vertex as a bad one */
  1024.             FLAG_VERTEX_AS_UNDEFINED(vlist[vert_free++]);
  1025.             continue;
  1026.         }
  1027.         map3d_xyz(points[i].x, points[i].y, points[i].z, vlist + vert_free);
  1028.         vlist[vert_free++].style = above;
  1029.         if (this_plot->plot_style == IMPULSES
  1030.             || this_plot->plot_style == BOXES) {
  1031.             map3d_xyz(points[i].x, points[i].y, min_array[FIRST_Z_AXIS], vlist + vert_free);
  1032.             vlist[vert_free++].style = above;
  1033.         }
  1034. #if FOUR_TRIANGLES
  1035.         /* FIXME: handle other styles as well! */
  1036.         if (this_plot->plot_style == LINES &&
  1037.             i < icrvs->p_count - 1)
  1038.             vert_free++;    /* keep one entry free for quad-center */
  1039. #endif
  1040.         }
  1041.     }
  1042.     }
  1043.  
  1044.     /* initialize plist: */
  1045.     id = 0;
  1046.     for (pfree = vert_free = 0, this_plot = plots, surface = 0;
  1047.      surface < pcount;
  1048.      this_plot = this_plot->next_sp, surface++) {
  1049.     row_offset = this_plot->iso_crvs->p_count;
  1050.     lp = &(this_plot->lp_properties);
  1051.     above = this_plot->lp_properties.l_type;
  1052.     below = this_plot->lp_properties.l_type + hiddenBacksideLinetypeOffset;
  1053.     GETNCRV(ncrv1);
  1054.     for (ncrv = 0, icrvs = this_plot->iso_crvs;
  1055.          ncrv < ncrv1 && icrvs;
  1056.          ncrv++, icrvs = icrvs->next)
  1057.         for (i = 0; i < icrvs->p_count; i++)
  1058.         switch (this_plot->plot_style) {
  1059.         case LINESPOINTS:
  1060.         case STEPS:    /* handle as LINES */
  1061.         case FSTEPS:
  1062.         case HISTEPS:
  1063.         case LINES:
  1064.             if (i < icrvs->p_count - 1 && ncrv < ncrv1 - 1) {
  1065. #if FOUR_TRIANGLES
  1066.             long *v3, *v4;
  1067.             int n3, n4;
  1068.             long line3, line4;
  1069.             int c3, c4, style3, style4;
  1070.             TBOOLEAN border3, border4, inc_id = FALSE;
  1071.  
  1072.             /* Build medium vertex: */
  1073.             vlist[vert_free + 1].x =
  1074.                 (vlist[vert_free].x + vlist[vert_free + 2 * row_offset - 1].x +
  1075.                  vlist[vert_free + 2].x + vlist[vert_free + 2 * row_offset + 1].x
  1076.                 ) / 4;
  1077.             vlist[vert_free + 1].y =
  1078.                 (vlist[vert_free].y + vlist[vert_free + 2 * row_offset - 1].y +
  1079.                  vlist[vert_free + 2].y + vlist[vert_free + 2 * row_offset + 1].y
  1080.                 ) / 4;
  1081.             vlist[vert_free + 1].z =
  1082.                 (vlist[vert_free].z + vlist[vert_free + 2 * row_offset - 1].z +
  1083.                  vlist[vert_free + 2].z + vlist[vert_free + 2 * row_offset + 1].z
  1084.                 ) / 4;
  1085.             vlist[vert_free + 1].style = above;
  1086.  
  1087.             PREPARE_POLYGON(n1, v1, 2, 0, 1, line1, c1,
  1088.                     border1, 0, -1, style1);
  1089.             /* !!FIXME!! special casing (see the older code) still needs to be
  1090.              * implemented! */
  1091.             if ((c1 || border1)
  1092.                 && reduce_polygon(&n1, &v1, &line1, 3)) {
  1093.                 CHECK_PLIST();
  1094.                 pfree += maybe_build_polygon(plist + pfree, n1, v1, line1,
  1095.                           style1, lp, -1L, -1L, id++,
  1096.                              is_untested);
  1097.                 inc_id = TRUE;
  1098.             }
  1099.             PREPARE_POLYGON(n2, v2, 0, 2 * row_offset - 1, 1, line2, c2,
  1100.                     border2, -1, 0, style2);
  1101.             if ((c2 || border2)
  1102.                 && reduce_polygon(&n2, &v2, &line2, 3)) {
  1103.                 CHECK_PLIST();
  1104.                 pfree += maybe_build_polygon(plist + pfree, n2, v2, line2,
  1105.                           style2, lp, -1L, -1L, id++,
  1106.                              is_untested);
  1107.                 inc_id = TRUE;
  1108.             }
  1109.             PREPARE_POLYGON(n3, v3, 2 * row_offset - 1, 2 * row_offset + 1, 1,
  1110.                   line3, c3, border3, icrvs->p_count - 2, -1,
  1111.                     style3);
  1112.             if ((c3 || border3)
  1113.                 && reduce_polygon(&n3, &v3, &line3, 3)) {
  1114.                 CHECK_PLIST();
  1115.                 pfree += maybe_build_polygon(plist + pfree, n3, v3, line3,
  1116.                           style3, lp, -1L, -1L, id++,
  1117.                              is_untested);
  1118.                 inc_id = TRUE;
  1119.             }
  1120.             PREPARE_POLYGON(n4, v4, 2 * row_offset + 1, 0, 1, line4, c4,
  1121.                     border4, -1, ncrv1 - 2, style4);
  1122.             if ((c4 || border4)
  1123.                 && reduce_polygon(&n4, &v4, &line4, 3)) {
  1124.                 CHECK_PLIST();
  1125.                 pfree += maybe_build_polygon(plist + pfree, n4, v4, line4,
  1126.                           style4, lp, -1L, -1L, id++,
  1127.                              is_untested);
  1128.                 inc_id = TRUE;
  1129.             }
  1130.             /*if (inc_id)
  1131.                id++; */
  1132.             vert_free++;
  1133.  
  1134. #else /* FOUR_TRIANGLES */
  1135.  
  1136.             PREPARE_POLYGON(n1, v1, row_offset, 0, 1, line1, c1,
  1137.                     border1, 0, 0, style1);
  1138.             PREPARE_POLYGON(n2, v2, 1, row_offset + 1, row_offset, line2,
  1139.                   c2, border2, icrvs->p_count - 2, ncrv1 - 2,
  1140.                     style2);
  1141.  
  1142.             /* if this makes at least one of the two triangles
  1143.              * visible, use the other diagonal here */
  1144.             if (hiddenShowAlternativeDiagonal
  1145.                 && !(v1)    /* HBB 980213: only try this in the case of */
  1146.                 &&!(v2)    /* undefined vertices -> *missing* trigs */
  1147.                 ) {
  1148.                 if (VERTEX_IS_UNDEFINED(vlist[vert_free + 1])) {
  1149.                 PREPARE_POLYGON(n1, v1, row_offset + 1, row_offset, 0,
  1150.                     line1, c1, border1, 0, ncrv1 - 2,
  1151.                         style1);
  1152.                 } else if (VERTEX_IS_UNDEFINED(vlist[vert_free + row_offset])) {
  1153.                 PREPARE_POLYGON(n1, v1, 0, 1, row_offset + 1, line1,
  1154.                       c1, border1, icrvs->p_count - 2, 0,
  1155.                         style1);
  1156.                 }
  1157.             }
  1158.             /* HBB 970323: if visible sides of the two triangles
  1159.              * differs, make diagonal visible! */
  1160.             /* HBB 970428: FIXME: should this be changed to only
  1161.              * activate *one* of the triangles' diagonals? If so, how
  1162.              * to find out the right one? Maybe the one with its
  1163.              * normal closer to the z direction? */
  1164.             if (hiddenHandleBentoverQuadrangles
  1165.                 && (c1 * c2) < 0
  1166.                 ) {
  1167.                 SETBIT(line1, n1 - 1, 1);
  1168.                 SETBIT(line2, n2 - 1, 1);
  1169.             }
  1170.             if ((c1 || border1)
  1171.                 && reduce_polygon(&n1, &v1, &line1, 3)) {
  1172.                 if ((c2 || border2)
  1173.                 && reduce_polygon(&n2, &v2, &line2, 3)) {
  1174.                 /* These two will need special handling, to
  1175.                  * ensure the links from one to the other
  1176.                  * are only set up when *both* polygons are
  1177.                  * valid */
  1178.                 int r1, r2;    /* temp. results */
  1179.  
  1180.                 CHECK_PLIST_2();
  1181.                 r1 = maybe_build_polygon(plist + pfree, n1, v1,
  1182.                           line1, style1, lp, -1L,
  1183.                            -1L, id, is_untested);
  1184.                 r2 = maybe_build_polygon(plist + pfree + r1, n2, v2,
  1185.                           line2, style2, lp, -1L,
  1186.                          -1L, id++, is_untested);
  1187.                 if (r1 && r2) {
  1188.                     plist[pfree].next_frag = pfree + 1;
  1189.                     plist[pfree + 1].next_frag = pfree;
  1190.                 }
  1191.                 pfree += r1 + r2;
  1192.                 } else {    /* P1 was ok, but P2 wasn't */
  1193. /* HBB 980213: P2 is invisible, so free its vertex list: */
  1194.                 if (v2) {
  1195.                     free(v2);
  1196.                     v2 = 0;
  1197.                 }
  1198. /* HBB 970323: if other polygon wasn't drawn, draw this polygon's
  1199.  * diagonal visibly (not only if it was 'on edge'). */
  1200.                 /*if (c2 || border2) */
  1201.                 SETBIT(line1, n1 - 1, 1);
  1202.                 CHECK_PLIST();
  1203.                 pfree += maybe_build_polygon(plist + pfree, n1, v1, line1,
  1204.                           style1, lp, -1L, -1L, id++,
  1205.                                  is_untested);
  1206.                 }
  1207.             } else {    /* P1 was not ok */
  1208. /* HBB 980213: P1 is invisible, so free its vertex list: */
  1209.                 if (v1) {
  1210.                 free(v1);
  1211.                 v1 = 0;
  1212.                 }
  1213. /* HBB 970323: if other polygon wasn't drawn, draw this polygon's
  1214.  * diagonal visibly (not only if it was 'on edge'). */
  1215.                 /*if (c1 || border1) */
  1216.                 SETBIT(line2, n2 - 1, 1);
  1217. /* HBB 970522: yet another problem triggered by the handling of
  1218.  * undefined points: if !(c2||border2), we mustn't try to reduce
  1219.  * triangle 2. */
  1220.                 if (0
  1221.                 || (1
  1222.                     && (c2 || border2)
  1223.                   && reduce_polygon(&n2, &v2, &line2, 2))
  1224.                 || (1
  1225. /* HBB 980213: reduce_... above might have zeroed v... */
  1226.                     && (v2)
  1227.                     && (c1 || border1))) {
  1228.                 CHECK_PLIST();
  1229.                 pfree += maybe_build_polygon(plist + pfree, n2, v2, line2,
  1230.                           style2, lp, -1L, -1L, id++,
  1231.                                  is_untested);
  1232.                 } else {
  1233.                 /* HBB 980213: sigh... both P1 and P2 were invisible, but
  1234.                  * at least one of them was a not undefined.. */
  1235.                 if (v1)
  1236.                     free(v1);
  1237.                 if (v2)
  1238.                     free(v2);
  1239.                 }
  1240.             }
  1241. #endif /* else: FOUR_TRIANGLES */
  1242.             }
  1243.             vert_free++;
  1244.             break;
  1245.         case DOTS:
  1246.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * 1,
  1247.                        "hidden v1 for dots");
  1248.             v1[0] = vert_free++;
  1249.             CHECK_PLIST();
  1250.             pfree += maybe_build_polygon(plist + pfree, 1, v1, 1L, above,
  1251.                     lp, -1L, -1L, id++, is_untested);
  1252.             break;
  1253.         case BOXERROR:    /* handle as POINTSTYLE */
  1254.         case BOXXYERROR:
  1255.         case XERRORBARS:    /* handle as POINTSTYLE */
  1256.         case YERRORBARS:
  1257.         case XYERRORBARS:
  1258.         case CANDLESTICKS:    /* HBB: these as well */
  1259.         case FINANCEBARS:
  1260.         case POINTSTYLE:
  1261.         case VECTOR:
  1262.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * 1,
  1263.                        "hidden v1 for point");
  1264.             v1[0] = vert_free++;
  1265.             CHECK_PLIST();
  1266.             pfree += maybe_build_polygon(plist + pfree, 1, v1, 0L, above,
  1267.                     lp, -1L, -1L, id++, is_untested);
  1268.             break;
  1269.         case BOXES:    /* handle as IMPULSES */
  1270.         case IMPULSES:
  1271.             n1 = 2;
  1272.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * n1,
  1273.                        "hidden v1 for impulse");
  1274.             v1[0] = vert_free++;
  1275.             v1[1] = vert_free++;
  1276.             line1 = 2L;
  1277.             (void) reduce_polygon(&n1, &v1, &line1, 1);
  1278.             CHECK_PLIST();
  1279.             pfree += maybe_build_polygon(plist + pfree, n1, v1, line1, above,
  1280.                     lp, -1L, -1L, id++, is_untested);
  1281.             break;
  1282.         }
  1283.     }
  1284. }
  1285.  
  1286. static int compare_by_zmax(p1, p2)
  1287. const void *p1, *p2;
  1288. {
  1289.     return (SIGNOF(plist[*(const long *) p2].zmax - plist[*(const long *) p1].zmax));
  1290. }
  1291.  
  1292. static void sort_by_zmax()
  1293.      /* and build the list (return_value = start of list) */
  1294. {
  1295.     long *sortarray, i;
  1296.     struct Polygon GPHUGE *this;
  1297.     sortarray = (long *) gp_alloc((unsigned long) sizeof(long) * pfree, "hidden sortarray");
  1298.     for (i = 0; i < pfree; i++)
  1299.     sortarray[i] = i;
  1300.     qsort(sortarray, (size_t) pfree, sizeof(long), compare_by_zmax);
  1301.     this = plist + sortarray[0];
  1302.     for (i = 1; i < pfree; i++) {
  1303.     this->next = sortarray[i];
  1304.     this = plist + sortarray[i];
  1305.     }
  1306.     this->next = -1L;
  1307.     pfirst = sortarray[0];
  1308.     free(sortarray);
  1309. }
  1310.  
  1311.  
  1312. /* HBB: try if converting this code to unsigned int (from signed shorts)
  1313.  *  fixes any of the remaining bugs. In the same motion, remove
  1314.  *  hardcoded sizeof(short) = 2 (09.10.1996) */
  1315.  
  1316. #define LASTBITS (USHRT_BITS -1)    /* ????? */
  1317. static int obscured(p)
  1318.      /* is p obscured by already drawn polygons? (only a simple minimax-test) */
  1319. struct Polygon GPHUGE *p;
  1320. {
  1321.     int l_xmin, l_xmax, l_ymin, l_ymax;        /* HBB 961110: avoid shadowing external names */
  1322.     t_pnt mask1, mask2;
  1323.     long indx1, indx2, k, m;
  1324.     tp_pnt cpnt;
  1325.     /* build the minimax-box */
  1326.     l_xmin = (p->xmin * xscaler) + xmiddle;
  1327.     l_xmax = (p->xmax * xscaler) + xmiddle;
  1328.     l_ymin = (p->ymin * yscaler) + ymiddle;
  1329.     l_ymax = (p->ymax * yscaler) + ymiddle;
  1330.     if (l_xmin < xleft)
  1331.     l_xmin = xleft;
  1332.     if (l_xmax > xright)
  1333.     l_xmax = xright;
  1334.     if (l_ymin < ybot)
  1335.     l_ymin = ybot;
  1336.     if (l_ymax > ytop)
  1337.     l_ymax = ytop;
  1338.     if (l_xmin > l_xmax || l_ymin > l_ymax)
  1339.     return 1;        /* not viewable */
  1340.     l_ymin = YREDUCE(l_ymin) - YREDUCE(ybot);
  1341.     l_ymax = YREDUCE(l_ymax) - YREDUCE(ybot);
  1342.     l_xmin = XREDUCE(l_xmin) - XREDUCE(xleft);
  1343.     l_xmax = XREDUCE(l_xmax) - XREDUCE(xleft);
  1344.     /* Now check bitmap */
  1345.     indx1 = l_ymin / PNT_BITS;
  1346.     indx2 = l_ymax / PNT_BITS;
  1347.     mask1 = PNT_MAX << (((unsigned) l_ymin) % PNT_BITS);
  1348.     mask2 = PNT_MAX >> ((~((unsigned) l_ymax)) % PNT_BITS);
  1349.     /* HBB: lclint wanted this: */
  1350.     assert(pnt != 0);
  1351.     for (m = l_xmin; m <= l_xmax; m++) {
  1352.     if (pnt[m] == 0)
  1353.         return 0;        /* not obscured */
  1354.     cpnt = pnt[m] + indx1;
  1355.     if (indx1 == indx2) {
  1356.         if (~(*cpnt) & mask1 & mask2)
  1357.         return 0;
  1358.     } else {
  1359.         if (~(*cpnt++) & mask1)
  1360.         return 0;
  1361.         for (k = indx1 + 1; k < indx2; k++)
  1362.         if ((*cpnt++) != PNT_MAX)
  1363.             return 0;
  1364.         if (~(*cpnt++) & mask2)
  1365.         return 0;
  1366.     }
  1367.     }
  1368.     return 1;
  1369. }
  1370.  
  1371.  
  1372. void draw_line_hidden(x1, y1, x2, y2)
  1373. unsigned int x1, y1, x2, y2;
  1374. {
  1375.     register struct termentry *t = term;
  1376.     TBOOLEAN flag;
  1377.     register int xv, yv, errx, erry, err;
  1378.     register unsigned int xvr, yvr;
  1379.     int unsigned xve, yve;
  1380.     register int dy, nstep = 0, dyr;
  1381.     int i;
  1382.  
  1383.     if (x1 > x2) {
  1384.     xvr = x2;
  1385.     yvr = y2;
  1386.     xve = x1;
  1387.     yve = y1;
  1388.     } else {
  1389.     xvr = x1;
  1390.     yvr = y1;
  1391.     xve = x2;
  1392.     yve = y2;
  1393.     };
  1394.     errx = XREDUCE(xve) - XREDUCE(xvr);
  1395.     erry = YREDUCE(yve) - YREDUCE(yvr);
  1396.     dy = (erry > 0 ? 1 : -1);
  1397.     dyr = dy * yfact;
  1398.     switch (dy) {
  1399.     case 1:
  1400.     nstep = errx + erry;
  1401.     errx = -errx;
  1402.     break;
  1403.     case -1:
  1404.     nstep = errx - erry;
  1405.     errx = -errx;
  1406.     erry = -erry;
  1407.     break;
  1408.     };
  1409.     err = errx + erry;
  1410.     errx <<= 1;
  1411.     erry <<= 1;
  1412.     xv = XREDUCE(xvr) - XREDUCE(xleft);
  1413.     yv = YREDUCE(yvr) - YREDUCE(ybot);
  1414.     (*t->move) (xvr, yvr);
  1415.     flag = !IS_UNSET(xv, yv);
  1416.     if (!hidden_no_update) {    /* Check first point */
  1417.     assert(ymax_hl != 0);
  1418.     assert(ymin_hl != 0);
  1419.     if (xv < xmin_hl)
  1420.         xmin_hl = xv;
  1421.     if (xv > xmax_hl)
  1422.         xmax_hl = xv;
  1423.     if (yv > ymax_hl[xv])
  1424.         ymax_hl[xv] = yv;
  1425.     if (yv < ymin_hl[xv])
  1426.         ymin_hl[xv] = yv;
  1427.     };
  1428.     for (i = 0; i < nstep; i++) {
  1429.     if (err < 0) {
  1430.         xv++;
  1431.         xvr += xfact;
  1432.         err += erry;
  1433.     } else {
  1434.         yv += dy;
  1435.         yvr += dyr;
  1436.         err += errx;
  1437.     };
  1438.     if (IS_UNSET(xv, yv)) {
  1439.         if (flag) {
  1440.         (*t->move) (xvr, yvr);
  1441.         flag = FALSE;
  1442.         };
  1443.     } else {
  1444.         if (!flag) {
  1445.         (*t->vector) (xvr, yvr);
  1446.         flag = TRUE;
  1447.         };
  1448.     };
  1449.     if (!hidden_no_update) {
  1450.         /* HBB 961110: lclint wanted these: */
  1451.         assert(ymax_hl != 0);
  1452.         assert(ymin_hl != 0);
  1453.         if (xv < xmin_hl)
  1454.         xmin_hl = xv;
  1455.         if (xv > xmax_hl)
  1456.         xmax_hl = xv;
  1457.         if (yv > ymax_hl[xv])
  1458.         ymax_hl[xv] = yv;
  1459.         if (yv < ymin_hl[xv])
  1460.         ymin_hl[xv] = yv;
  1461.     };
  1462.     };
  1463.     if (!flag)
  1464.     (*t->vector) (xve, yve);
  1465.     return;
  1466. }
  1467.  
  1468.  
  1469. static GP_INLINE int xy_overlap(a, b)
  1470.      /* Do a and b overlap in x or y (minimax test) */
  1471. struct Polygon GPHUGE *a, GPHUGE * b;
  1472. {
  1473. #if TEST_GRIDCHECK
  1474.     /* First, check by comparing the extent bit patterns: */
  1475.     if (!((a->xextent & b->xextent) && (a->yextent & b->yextent)))
  1476.     return 0;
  1477. #endif
  1478.     return ((int) (GE(b->xmax, a->xmin) && GE(a->xmax, b->xmin)
  1479.            && GE(b->ymax, a->ymin) && GE(a->ymax, b->ymin)));
  1480. }
  1481.  
  1482.  
  1483. static void get_plane(p, a, b, c, d)
  1484. struct Polygon GPHUGE *p;
  1485. double *a, *b, *c, *d;
  1486. {
  1487.     int i;
  1488.     struct Vertex GPHUGE *v1, GPHUGE * v2;
  1489.     double x, y, z, s;
  1490.     if (p->n == 1) {
  1491.     *a = 0.0;
  1492.     *b = 0.0;
  1493.     *c = 1.0;
  1494.     *d = -vlist[p->vertex[0]].z;
  1495.     return;
  1496.     }
  1497.     v1 = vlist + p->vertex[p->n - 1];
  1498.     v2 = vlist + p->vertex[0];
  1499.     *a = (v1->y - v2->y) * (v1->z + v2->z);
  1500.     *b = (v1->z - v2->z) * (v1->x + v2->x);
  1501.     *c = (v1->x - v2->x) * (v1->y + v2->y);
  1502.     for (i = 0; i < p->n - 1; i++) {
  1503.     v1 = vlist + p->vertex[i];
  1504.     v2 = vlist + p->vertex[i + 1];
  1505.     *a += (v1->y - v2->y) * (v1->z + v2->z);
  1506.     *b += (v1->z - v2->z) * (v1->x + v2->x);
  1507.     *c += (v1->x - v2->x) * (v1->y + v2->y);
  1508.     }
  1509.     s = sqrt(*a * *a + *b * *b + *c * *c);
  1510.     if (GE(0.0, s)) {        /* => s==0 => the vertices are in one line */
  1511.     v1 = vlist + p->vertex[0];
  1512.     for (i = 1; i < p->n; i++) {
  1513.         v2 = vlist + p->vertex[i];
  1514.         if (!V_EQUAL(v1, v2))
  1515.         break;
  1516.     }
  1517.     /* (x,y,z) is l.u. from <v1, v2> */
  1518.     x = v1->x;
  1519.     y = v1->y;
  1520.     z = v1->z;
  1521.     if (EQ(y, v2->y))
  1522.         y += 1.0;
  1523.     else
  1524.         x += 1.0;
  1525.     /* build a vector that is orthogonal to the line of the polygon */
  1526.     *a = v1->y * (v2->z - z) + v2->y * (z - v1->z) + y * (v1->z - v2->z);
  1527.     *b = v1->z * (v2->x - x) + v2->z * (x - v1->x) + z * (v1->x - v2->x);
  1528.     *c = v1->x * (v2->y - y) + v2->x * (y - v1->y) + x * (v1->y - v2->y);
  1529.     s = sqrt(*a * *a + *b * *b + *c * *c);
  1530.     }
  1531.     if (*c < 0.0)        /* ensure that normalized c is > 0 */
  1532.     s *= -1.0;        /* The exact relation, because c can be very small */
  1533.     *a /= s;
  1534.     *b /= s;
  1535.     *c /= s;
  1536.     *d = -*a * v1->x - *b * v1->y - *c * v1->z;
  1537.     return;
  1538. }
  1539.  
  1540.  
  1541. static t_poly_plane_intersect
  1542.  polygon_plane_intersection(p, a, b, c, d)
  1543. struct Polygon GPHUGE *p;
  1544. double a, b, c, d;
  1545. {
  1546.     int i, sign, max_sign, min_sign;
  1547.     struct Vertex GPHUGE *v;
  1548.  
  1549.     CHECK_POINTER(plist, p);
  1550.  
  1551.     v = vlist + p->vertex[0];
  1552.     max_sign = min_sign = SIGNOF(a * v->x + b * v->y + c * v->z + d);
  1553.     for (i = 1; i < p->n; i++) {
  1554.     v = vlist + p->vertex[i];
  1555.     sign = SIGNOF(a * v->x + b * v->y + c * v->z + d);
  1556.     if (sign > max_sign)
  1557.         max_sign = sign;
  1558.     else if (sign < min_sign)
  1559.         min_sign = sign;
  1560.     }
  1561.  
  1562.     /* Yes, this is a bit tricky, but it's simpler for the computer... */
  1563.     if (min_sign == -1) {
  1564.     if (max_sign == 1)
  1565.         return (intersects);
  1566.     else
  1567.         return behind;
  1568.     } else {
  1569.     if ((max_sign == 0) && (min_sign == 0))
  1570.         return (inside);
  1571.     else
  1572.         return infront;
  1573.     }
  1574. }
  1575.  
  1576.  
  1577. /* full xy-overlap test (support routines first) */
  1578.  
  1579. /* What do negative return values mean?
  1580.    It is allowed, that the 2 polygons touch themselves in one point.
  1581.    There are 2 possibilities of this case:
  1582.    (A) A vertex of the one polygon is on an edge of the other
  1583.    (B) The two polygons have a vertex together.
  1584.    In case (A) the algorithm finds 2 edges, wich touches an edge of the
  1585.    other polygon. In case (B) the algorithm finds four pairs of edges.
  1586.    I handle both cases with negative return values:
  1587.    Case (A) returns -2, and case (B) returns -1.
  1588.    So I can say, if the sum of negative return values goes greater than 4
  1589.    (absolutly), there must be more than one touch point. So the 
  1590.    polygons must overlap. */
  1591.  
  1592. /* some variables, for keeping track of minima and maxima across
  1593.  * all these routines */
  1594.  
  1595. static double m1x, M1x, m1y, M1y, m2x, M2x, m2y, M2y;
  1596.  
  1597. #define MINIMAX  (GE(M2x,m1x) && GE(M1x,m2x) && GE(M2y,m1y) && GE(M1y,m2y))
  1598.  
  1599.  
  1600. /* Does the edge [v1,v2] intersect edge [v3, v4] ?
  1601.  * This is for non-vertical [v1,v2]
  1602.  */
  1603. static int intersect(v1, v2, v3, v4)
  1604. struct Vertex GPHUGE *v1, GPHUGE * v2, GPHUGE * v3, GPHUGE * v4;
  1605. {
  1606.     double m1, m2, t1, t2, x, y, minx, maxx;
  1607.  
  1608.     m1 = (v2->y - v1->y) / (v2->x - v1->x);
  1609.     t1 = v1->y - m1 * v1->x;
  1610.  
  1611.     /* if [v3,v4] vertical: */
  1612.     if (EQ(v3->x, v4->x)) {
  1613.     y = v3->x * m1 + t1;
  1614.     if (GR(m2x, m1x) && GR(M1x, m2x)) {
  1615.         if (GR(y, m2y) && GR(M2y, y))
  1616.         return 1;
  1617.         if (EQ(y, m2y) || EQ(y, M2y))
  1618.         return -2;
  1619.         return 0;
  1620.     }
  1621.     /* m2x == m1x || m2x == M1x */
  1622.     if (GR(y, m2y) && GR(M2y, y))
  1623.         return -2;
  1624.     if (EQ(y, m2y) || EQ(y, M2y))
  1625.         return -1;
  1626.     return 0;
  1627.     }
  1628.     /* [v3,v4] not vertical */
  1629.     m2 = (v4->y - v3->y) / (v4->x - v3->x);
  1630.     t2 = v3->y - m2 * v3->x;
  1631.     if (!SIGNOF(m1 - m2)) {    /* if edges are parallel: */
  1632.     x = m1 * v3->x + t1 - v3->y;
  1633.     if (!EQ(m1 * v3->x + t1 - v3->y, 0.0))
  1634.         return 0;
  1635.     if (GR(M2x, m1x) && GR(M1x, m2x))
  1636.         return 1;
  1637.     return -1;        /* the edges have a common vertex */
  1638.     }
  1639.     x = (t2 - t1) / (m1 - m2);
  1640.     minx = GPMAX(m1x, m2x);
  1641.     maxx = GPMIN(M1x, M2x);
  1642.     if (GR(x, minx) && GR(maxx, x))
  1643.     return 1;
  1644.     if (GR(minx, x) || GR(x, maxx))
  1645.     return 0;
  1646.     if ((EQ(x, m1x) || EQ(x, M1x)) && (EQ(x, m1x) || EQ(x, M1x)))
  1647.     return -1;
  1648.     return -2;
  1649. }
  1650.  
  1651. /* Does the edge [v1,v2] intersect edge [v3, v4] ?
  1652.  * This is for vertical [v1,v2]
  1653.  */
  1654. static int v_intersect(v1, v2, v3, v4)
  1655. struct Vertex GPHUGE *v1, GPHUGE * v2, GPHUGE * v3, GPHUGE * v4;
  1656. {
  1657.     double y;
  1658.  
  1659.     /* HBB 971115: to silence picky compilers, use v2 at least once, in
  1660.      * a dummy expression (caution, may be ANSI-C specific because of
  1661.      * the use of 'void'...) */
  1662.     (void) v2;
  1663.  
  1664.     /* if [v3,v4] is vertical, too: */
  1665.     /* already known: rectangles do overlap, because MINIMAX was true */
  1666.     if (EQ(v3->x, v4->x))
  1667.     return (GR(M2y, m1y) && GR(M1y, m2y)) ? 1 : -1;
  1668.  
  1669.     /*  [v3,v4] not vertical */
  1670.     y = v3->y + (v1->x - v3->x) * (v4->y - v3->y) / (v4->x - v3->x);
  1671.     if (GR(m1x, m2x) && GR(M2x, m1x)) {
  1672.     if (GR(y, m1y) && GR(M1y, y))
  1673.         return 1;
  1674.     if (EQ(y, m1y) || EQ(y, M1y))
  1675.         return -2;
  1676.     return 0;
  1677.     }
  1678.     /* m1x == m2x || m1x == M2x */
  1679.     if (GR(y, m1y) && GR(M1y, y))
  1680.     return -2;
  1681.     if (EQ(y, m1y) || EQ(y, M1y))
  1682.     return -1;
  1683.     return 0;
  1684. }
  1685.  
  1686. #define UPD_MINMAX(v1,v2,npar) do {    \
  1687.   if (v1->x< v2->x)                 \
  1688.     CONCAT3(m,npar,x) = v1->x, CONCAT3(M,npar,x) = v2->x;   \
  1689.   else                              \
  1690.     CONCAT3(m,npar,x) = v2->x, CONCAT3(M,npar,x) = v1->x;   \
  1691.   if (v1->y< v2->y)                 \
  1692.     CONCAT3(m,npar,y) = v1->y, CONCAT3(M,npar,y) = v2->y;   \
  1693.   else                              \
  1694.     CONCAT3(m,npar,y) = v2->y, CONCAT3(M,npar,y) = v1->y;   \
  1695. } while (0)
  1696.  
  1697. /* does the edge [v1,v2] intersect polygon p? */
  1698. static int intersect_polygon(v1, v2, p)
  1699. struct Vertex GPHUGE *v1, GPHUGE * v2;
  1700. struct Polygon GPHUGE *p;
  1701. {
  1702.     struct Vertex GPHUGE *v3, GPHUGE * v4;
  1703.     int i, s, t = 0;
  1704.     int (*which_intersect) __PROTO((struct Vertex GPHUGE * v1,
  1705.                     struct Vertex GPHUGE * v2,
  1706.                     struct Vertex GPHUGE * v3,
  1707.                     struct Vertex GPHUGE * v4))
  1708.     = intersect;
  1709.  
  1710.     /* Is [v1,v2] vertical? If, use v_intersect() */
  1711.     if (EQ(v1->x, v2->x))
  1712.     which_intersect = v_intersect;
  1713.  
  1714.     UPD_MINMAX(v1,v2,1);
  1715.  
  1716.     /* test first edge of polygon p */
  1717.     v3 = vlist + p->vertex[p->n - 1];
  1718.     v4 = vlist + p->vertex[0];
  1719.     UPD_MINMAX(v3,v4,2);
  1720.  
  1721.  
  1722.     if (MINIMAX) {
  1723.     s = (*which_intersect) (v1, v2, v3, v4);
  1724.     if (s == 1 || (s < 0 && (t += s) < -4))
  1725.         return 1;
  1726.     }
  1727.     /* and the other edges... */
  1728.     for (i = 0; i < p->n - 1; i++) {
  1729.     v3 = vlist + p->vertex[i];
  1730.     v4 = vlist + p->vertex[i + 1];
  1731.     UPD_MINMAX(v3,v4,2);
  1732.     if (!MINIMAX)
  1733.         continue;
  1734.     s = (*which_intersect) (v1, v2, v3, v4);
  1735.     if (s == 1 || (s < 0 && (t += s) < -4))
  1736.         return 1;
  1737.     }
  1738.     return t;
  1739. }
  1740.  
  1741. /* OK, now here comes the 'real' polygon intersection routine:
  1742.  * Do a and b overlap in x or y (full test):
  1743.  * Do edges intersect, do the polygons touch in two points,
  1744.  * or is a vertex of one polygon inside the other polygon? */
  1745.  
  1746. static int full_xy_overlap(a, b)
  1747. struct Polygon GPHUGE *a, GPHUGE * b;
  1748. {
  1749.     struct Vertex GPHUGE *v1, GPHUGE * v2, v;
  1750.     int i, s, t = 0;
  1751.  
  1752.     if (a->n < 2 || b->n < 2)
  1753.     return 0;
  1754.     v1 = vlist + a->vertex[a->n - 1];
  1755.     v2 = vlist + a->vertex[0];
  1756.     s = intersect_polygon(v1, v2, b);
  1757.     if (s == 1 || (s < 0 && (t += s) < -4))
  1758.     return 1;
  1759.     for (i = 0; i < a->n - 1; i++) {
  1760.     v1 = vlist + a->vertex[i];
  1761.     v2 = vlist + a->vertex[i + 1];
  1762.     s = intersect_polygon(v1, v2, b);
  1763.     if (s == 1 || (s < 0 && (t += s) < -4))
  1764.         return 1;
  1765.     }
  1766.     /* No edges intersect. Is one polygon inside the other? */
  1767.     /* The  inner polygon has the greater minx */
  1768.     if (a->xmin < b->xmin) {
  1769.     struct Polygon GPHUGE *temp = a;
  1770.     a = b;
  1771.     b = temp;
  1772.     }
  1773.     /* Now, a is the inner polygon */
  1774.     for (i = 0; i < a->n; i++) {
  1775.     v1 = vlist + a->vertex[i];
  1776.     /* HBB: lclint seems to have found a problem here: v wasn't
  1777.      * fully defined when passed to intersect_polygon() */
  1778.     v = *v1;
  1779.     v.x = -1.1;
  1780.     s = intersect_polygon(v1, &v, b);
  1781.     if (s == 0)
  1782.         return 0;        /* a is outside polygon b */
  1783.     v.x = 1.1;
  1784.     if (s == 1) {
  1785.         s = intersect_polygon(v1, &v, b);
  1786.         if (s == 0)
  1787.         return 0;
  1788.         if (s == 1)
  1789.         return 1;
  1790.     } else if (intersect_polygon(v1, &v, b) == 0)
  1791.         return 0;
  1792.     }
  1793. #if 1
  1794.     /* 970614: don't bother, just presume that they do overlap: */
  1795.     return 1;
  1796. #else /* !1 */
  1797.     print_polygon(a, "a");
  1798.     print_polygon(b, "b");
  1799.     graph_error("Logic Error in full_xy_overlap");
  1800.     return -1;            /* HBB: shut up gcc -Wall */
  1801. #endif /* !1 */
  1802. }
  1803.  
  1804.  
  1805. /* Helper routine for split_polygon_by_plane */
  1806. static long build_new_vertex(V1, V2, w)
  1807. long V1, V2;    /* the vertices, between which a new vertex is demanded */
  1808. double w;    /* information about where between V1 and V2 it should be */
  1809. {
  1810.     long V;
  1811.     struct Vertex GPHUGE *v, GPHUGE * v1, GPHUGE * v2;
  1812.  
  1813.     if (EQ(w, 0.0))
  1814.     return V1;
  1815.     if (EQ(w, 1.0))
  1816.     return V2;
  1817.  
  1818.     /* We need a new Vertex */
  1819.     if (vert_free == nvert)    /* Extend vlist, if necessary */
  1820.     vlist = (struct Vertex GPHUGE *)
  1821.         gp_realloc(vlist, (unsigned long) sizeof(struct Vertex) * (nvert += 10L), "hidden vlist");
  1822.     V = vert_free++;
  1823.     v = vlist + V;
  1824.     v1 = vlist + V1;
  1825.     v2 = vlist + V2;
  1826.  
  1827.     v->x = (v2->x - v1->x) * w + v1->x;
  1828.     v->y = (v2->y - v1->y) * w + v1->y;
  1829.     v->z = (v2->z - v1->z) * w + v1->z;
  1830.     v->style = -1;
  1831.  
  1832.     /* 970610: additional checks, to prevent adding unnecessary vertices */
  1833.     if (V_EQUAL(v, v1)) {
  1834.     vert_free--;
  1835.     return V1;
  1836.     }
  1837.     if (V_EQUAL(v, v2)) {
  1838.     vert_free--;
  1839.     return V2;
  1840.     }
  1841.     return V;
  1842. }
  1843.  
  1844.  
  1845. /* Splits polygon p by the plane represented by its equation
  1846.  * coeffecients a to d.
  1847.  * return-value: part of the polygon, that is in front/back of plane
  1848.  * (Distinction necessary to ensure the ordering of polygons in
  1849.  * the plist stays intact)
  1850.  * Caution: plist and vlist can change their location!!!
  1851.  * If a point is in plane, it comes to the behind - part
  1852.  * (HBB: that may well have changed, as I cut at the 'in plane'
  1853.  * mechanisms heavily)
  1854.  */
  1855.  
  1856. static long split_polygon_by_plane(P, a, b, c, d, f)
  1857. long P;                /* Polygon as index on plist */
  1858. double a, b, c, d;
  1859. TBOOLEAN f;            /* return value = Front(1) or Back(0) */
  1860. {
  1861.     int i, j;
  1862.     struct Polygon GPHUGE *p = plist + P;
  1863.     struct Vertex GPHUGE *v;
  1864.     int snew, stest;
  1865.     int in_plane;        /* what is about points in plane? */
  1866.     int cross1, cross2;        /* first vertices, after p crossed the plane */
  1867.     double a1 = 0.0, b1, a2 = 0.0, b2;    /* parameters of the crosses */
  1868.     long Front;            /* the new Polygon */
  1869.     int n1, n2;
  1870.     long *v1, *v2;
  1871.     long line1, line2;
  1872.     long next_frag_temp;
  1873.  
  1874.  
  1875.     CHECK_POINTER(plist, p);
  1876.  
  1877.     in_plane = (EQ(c, 0.0) && f) ? 1 : -1;
  1878.     /* first vertex */
  1879.     v = vlist + p->vertex[0];
  1880.  
  1881.     CHECK_POINTER(vlist, v);
  1882.  
  1883.     b1 = a * v->x + b * v->y + c * v->z + d;
  1884.     stest = SIGNOF(b1);
  1885. #if OLD_SPLIT_INPLANE
  1886.     if (!stest)
  1887.     stest = in_plane;
  1888. #endif
  1889.  
  1890.     /* find first vertex that is on the other side of the plane */
  1891.     for (cross1 = 1, snew = stest; snew == stest && cross1 < p->n; cross1++) {
  1892.     a1 = b1;
  1893.     v = vlist + p->vertex[cross1];
  1894.     CHECK_POINTER(vlist, v);
  1895.     b1 = a * v->x + b * v->y + c * v->z + d;
  1896.     snew = SIGNOF(b1);
  1897. #if OLD_SPLIT_INPLANE
  1898.     if (!snew)
  1899.         snew = in_plane;
  1900. #endif
  1901.     }
  1902.  
  1903.     if (snew == stest) {
  1904.     /*HBB: this is possibly not an error in split_polygon, it's just
  1905.      * not possible to split this polygon by this plane. I.e., the
  1906.      * error is in the caller!  */
  1907.     fprintf(stderr, "split_poly failed, polygon nr. %ld\n", P);
  1908.     return (-1);        /* return 'unsplittable' */
  1909.     }
  1910.     cross1--;            /* now, it's the last one on 'this' side */
  1911.  
  1912.     /* first vertex that is on the first side again */
  1913.     for (b2 = b1, cross2 = cross1 + 1;
  1914.      snew != stest && cross2 < p->n; cross2++) {
  1915.     a2 = b2;
  1916.     v = vlist + p->vertex[cross2];
  1917.     CHECK_POINTER(vlist, v);
  1918.     b2 = a * v->x + b * v->y + c * v->z + d;
  1919.     snew = SIGNOF(b2);
  1920. #if OLD_SPLIT_INPLANE
  1921.     if (!snew)
  1922.         snew = -in_plane;
  1923. #endif
  1924.     }
  1925.  
  1926.     if (snew != stest) {
  1927.     /* only p->vertex[0] is on 'this' side */
  1928.     a2 = b2;
  1929.     v = vlist + p->vertex[0];
  1930.     CHECK_POINTER(vlist, v);
  1931.     b2 = a * v->x + b * v->y + c * v->z + d;
  1932.     } else
  1933.     cross2--;        /* now it's the last one on 'the other' side */
  1934.  
  1935.     /* We need two new polygons instead of the old one: */
  1936.     n1 = p->n - cross2 + cross1 + 2;
  1937.     n2 = cross2 - cross1 + 2;
  1938.     v1 = (long *) gp_alloc((unsigned long) sizeof(long) * n1,
  1939.                "hidden v1 for two new poly");
  1940.     v2 = (long *) gp_alloc((unsigned long) sizeof(long) * n2,
  1941.                "hidden v2 for two new poly");
  1942. #if SHOW_SPLITTING_EDGES
  1943.     line1 = 1L << (n1 - 1);
  1944.     line2 = 1L << (n2 - 1);
  1945. #else
  1946.     line1 = line2 = 0L;
  1947. #endif
  1948.     v1[0] = v2[n2 - 1] =
  1949.     build_new_vertex(p->vertex[cross2 - 1],
  1950.         p->vertex[(cross2 < p->n) ? cross2 : 0], a2 / (a2 - b2));
  1951.     v2[0] = v1[n1 - 1] =
  1952.     build_new_vertex(p->vertex[cross1 - 1],
  1953.              p->vertex[cross1], a1 / (a1 - b1));
  1954.  
  1955.     /* Copy visibility from split edges to their fragments: */
  1956.     if (p->line & (1 << (cross1 - 1))) {
  1957.     line1 |= 1L << (n1 - 2);
  1958.     line2 |= 1L;
  1959.     }
  1960.     if (p->line & (1L << (cross2 - 1))) {
  1961.     line1 |= 1L;
  1962.     line2 |= 1L << (n2 - 2);
  1963.     }
  1964.     /* Copy the old line patterns and vertex listings into new places */
  1965.     /* p->vertex[cross2...p->n-1] --> v1[1...] */
  1966.     for (i = cross2, j = 1; i < p->n;) {
  1967.     if (p->line & (1L << i))
  1968.         line1 |= 1L << j;
  1969.     v1[j++] = p->vertex[i++];
  1970.     }
  1971.     /* p->vertex[0...cross1-1] --> v1[...n1-2] */
  1972.     for (i = 0; i < cross1 - 1;) {
  1973.     if (p->line & (1L << i))
  1974.         line1 |= 1L << j;
  1975.     v1[j++] = p->vertex[i++];
  1976.     }
  1977.     v1[j++] = p->vertex[i++];
  1978.  
  1979.     /* p->vertex[cross1...cross2-1] -> v2[1...n2-2] */
  1980.     for (j = 1; i < cross2 - 1;) {
  1981.     if (p->line & (1L << i))
  1982.         line2 |= 1L << j;
  1983.     v2[j++] = p->vertex[i++];
  1984.     }
  1985.     v2[j++] = p->vertex[i];
  1986.  
  1987.     /* old vertex list isn't needed any more: */
  1988.     free(p->vertex);
  1989.     p->vertex = 0;
  1990.  
  1991.     if (!reduce_polygon(&n1, &v1, &line1, 1) ||
  1992.     !reduce_polygon(&n2, &v2, &line2, 1))
  1993.     graph_error("Logic Error 2 in split_polygon");
  1994.  
  1995.     /* Build the 2 new polygons : we are reusing one + making one new */
  1996.     CHECK_PLIST_EXTRA(p = plist + P);
  1997.     Front = pfree++;
  1998.  
  1999.     if ((next_frag_temp = p->next_frag) < 0)
  2000.     next_frag_temp = P;    /* First split of this polygon at all */
  2001.  
  2002.     /* HBB 961110: lclint said == shouldn't be used on Boolean's */
  2003.     if ((f && (stest < 0)) || ((!f) && !(stest < 0))) {
  2004.     build_polygon(plist + P, n1, v1, line1, p->style, p->lp,
  2005.               p->next, Front, p->id, p->tested);
  2006.     build_polygon(plist + Front, n2, v2, line2, p->style, p->lp,
  2007.               -1, next_frag_temp, p->id, is_moved_or_split);
  2008.     } else {
  2009.     build_polygon(plist + P, n2, v2, line2, p->style, p->lp,
  2010.               p->next, Front, p->id, p->tested);
  2011.     build_polygon(plist + Front, n1, v1, line1, p->style, p->lp,
  2012.               -1, next_frag_temp, p->id, is_moved_or_split);
  2013.     }
  2014.     return Front;
  2015. }
  2016.  
  2017. /* HBB: these pieces of code are found repeatedly in in_front(), so
  2018.  * I put them into macros
  2019.  * Note: can't use the 'do {...} while (0)' method for
  2020.  * these: 'continue' wouldn't work any more.
  2021.  */
  2022. #define PUT_P_IN_FRONT_TEST(new_tested) {\
  2023.   plist[Plast].next = p->next; \
  2024.   p->next = Test; \
  2025.   p->tested = new_tested; \
  2026.   if (Insert >= 0) \
  2027.     plist[Insert].next = P; \
  2028.   else \
  2029.     pfirst = P; \
  2030.   Insert = P; \
  2031.   P = Plast; \
  2032.   p = plist + P; \
  2033.   continue; \
  2034. }
  2035.  
  2036. #define SPLIT_TEST_BY_P {\
  2037.   long Back = split_polygon_by_plane (Test, pa, pb, pc, pd, 0); \
  2038.   if (Back >0) { \
  2039.     p = plist + P; \
  2040.     test = plist + Test;  /* plist can change */ \
  2041.     plist[Back].next = p->next; \
  2042.     p->next = Back; \
  2043.     P = Back; \
  2044.     p = plist + P; \
  2045.     zmin = test->zmin;  /* the new z-value */ \
  2046.   } else { \
  2047.     fprintf(stderr, "Failed to split test (%ld) by p (%ld)\n", Test, P); \
  2048.     graph_error("Error in hidden line removal: couldn't split..."); \
  2049.   } \
  2050.   continue; \
  2051. }
  2052.  
  2053. #define SPLIT_P_BY_TEST {\
  2054.   long Front = split_polygon_by_plane (P, ta, tb, tc, td, 1);\
  2055.   if (Front >0) {\
  2056.     p = plist + P;\
  2057.     test = plist + Test;    /* plist can change */\
  2058.     plist[Front].next = Test;\
  2059.     if (Insert >= 0)\
  2060.       plist[Insert].next = Front;\
  2061.     else\
  2062.       pfirst = Front;\
  2063.     Insert = Front;\
  2064.   } else {\
  2065.     fprintf(stderr, "Failed to split p (%ld) by test(%ld), relations are %d, %d\n", P, Test, p_rel_tplane, t_rel_pplane);\
  2066.     print_polygon(test, "test");\
  2067.     print_polygon(p, "p");\
  2068.     fputc('\n', stderr);\
  2069.   }\
  2070.   continue; /* FIXME: should we continue nevertheless? */\
  2071. }
  2072.  
  2073.  
  2074. /* Is Test in front of the other polygons?
  2075.  * If not, polygons which are in front of test come between
  2076.  * Last and Test (I.e. to the 'front' of the plist)
  2077.  */
  2078. static int in_front(Last, Test)
  2079. long Last, Test;
  2080. {
  2081.     struct Polygon GPHUGE *test = plist + Test, GPHUGE * p;
  2082.     long Insert, P, Plast;
  2083.     double zmin,        /* minimum z-value of Test */
  2084.      ta, tb, tc, td,        /* the plane of Test      */
  2085.      pa, pb, pc, pd;        /* the plane of a polygon */
  2086.     t_poly_plane_intersect p_rel_tplane, t_rel_pplane;
  2087.     int loop = (test->tested == is_in_loop);
  2088.  
  2089.     CHECK_POINTER(plist, test);
  2090.  
  2091.     /* Maybe this should only done at the end of this routine... */
  2092.     /* test->tested = is_tested; */
  2093.  
  2094.     Insert = Last;
  2095.  
  2096.     /* minimal z-value of Test */
  2097.     zmin = test->zmin;
  2098.  
  2099.     /* The plane of the polygon Test */
  2100.     get_plane(test, &ta, &tb, &tc, &td);
  2101.  
  2102.     /* Compare Test with the following polygons, which overlap in z value */
  2103.     for (Plast = Test, p = test;
  2104.      ((P = p->next) >= 0) && (((p = plist + P)->zmax > zmin) || (p->tested != 0));
  2105.      Plast = P) {
  2106.     CHECK_POINTER(plist, p);
  2107.  
  2108.     if (!xy_overlap(test, p))
  2109.         continue;
  2110.  
  2111.     if (p->zmax <= zmin)
  2112.         continue;
  2113.  
  2114.     p_rel_tplane = polygon_plane_intersection(p, ta, tb, tc, td);
  2115.     if ((p_rel_tplane == behind) || (p_rel_tplane == inside))
  2116.         continue;
  2117.  
  2118.     get_plane(p, &pa, &pb, &pc, &pd);
  2119.     t_rel_pplane = polygon_plane_intersection(test, pa, pb, pc, pd);
  2120.     if ((t_rel_pplane == infront) || (t_rel_pplane == inside))
  2121.         continue;
  2122.  
  2123.     if (!full_xy_overlap(test, p))
  2124.         continue;
  2125.  
  2126. #if 1
  2127.     /* HBB 971115: try new approach developed directly from Foley et
  2128.      * al.'s original description */
  2129.     /* OK, at this point, a total of 16 cases is left to be handled,
  2130.      * as 4 variables are important, each with two possible values:
  2131.      * t_rel_pplane may be 'behind' or 'intersects', p_rel_tplane may
  2132.      * be 'infront' or 'intersects', 'test' may be (seem to be) part
  2133.      * of a loop or not ('loop'), and p may have been considered
  2134.      * already or not (p->tested =?= is_tested). */
  2135.     /* See if splitting wherever possible is a good idea: */
  2136.     if (p_rel_tplane == intersects) {
  2137.         p->tested = is_moved_or_split;
  2138.         SPLIT_P_BY_TEST;
  2139.     } else if (t_rel_pplane == intersects) {
  2140.         test->tested = is_moved_or_split;
  2141.         SPLIT_TEST_BY_P;
  2142.     } else {
  2143.         if (loop && (p->tested == is_tested)) {
  2144.         /* Ouch, seems like we're in trouble, really */
  2145.         fprintf(stderr, "\
  2146. #Failed: In loop, without intersections.\n\
  2147. #Relations are %d, %d\n",
  2148.             p_rel_tplane, t_rel_pplane);
  2149.         print_polygon(test, "test");
  2150.         print_polygon(p, "p");
  2151.         continue;    /* Keep quiet, maybe no-one will notice (;-) */
  2152.         } else {
  2153.         PUT_P_IN_FRONT_TEST(is_in_loop);
  2154.         }
  2155.     }
  2156. #else
  2157.     /* p obscures test (at least, this *might* be happening */
  2158.     if (loop) {
  2159.         /* if P isn't interesting for the loop, put P in front of Test */
  2160.         if ((p->tested != is_tested)
  2161.         || (p_rel_tplane == infront)
  2162.         || (t_rel_pplane == behind))    /* HBB 970325: was infront (!?) */
  2163.         PUT_P_IN_FRONT_TEST(is_moved_or_split);
  2164.  
  2165.         /* else, split either test or p */
  2166.         if (t_rel_pplane == intersects)
  2167.         SPLIT_TEST_BY_P;
  2168.         if (p_rel_tplane == intersects)
  2169.         SPLIT_P_BY_TEST;
  2170.  
  2171.         /* HBB: if we ever come through here, something went severely wrong */
  2172.         fprintf(stderr, "Failed: polygon %ld vs. %ld, relations are %d, %d\n", P, Test, p_rel_tplane, t_rel_pplane);
  2173.         print_polygon(test, "test");
  2174.         print_polygon(p, "p");
  2175.         fputc('\n', stderr);
  2176.         graph_error("Couldn't resolve a polygon overlapping pb. Go tell HBB...");
  2177.     }
  2178.     /* No loop: if it makes sense, put P in front of Test */
  2179.     if ((p->tested != is_tested)
  2180.         && ((t_rel_pplane == behind)
  2181.         || (p_rel_tplane == infront)))
  2182.         PUT_P_IN_FRONT_TEST(is_moved_or_split);
  2183.  
  2184.     /* If possible, split P */
  2185.     if ((p->tested != is_tested) || (p_rel_tplane == intersects))
  2186.         SPLIT_P_BY_TEST;
  2187.  
  2188.     /* if possible, split Test */
  2189.     if (t_rel_pplane == intersects)
  2190.         SPLIT_TEST_BY_P;
  2191.  
  2192.     /* else, we have a loop: mark P as loop and put it in front of Test */
  2193.     PUT_P_IN_FRONT_TEST(is_in_loop);    /* -2: p is in a loop */
  2194. #endif
  2195.     }
  2196.     if (Insert == Last)
  2197.     test->tested = is_tested;
  2198.     return (int) (Insert == Last);
  2199. }
  2200.  
  2201.  
  2202. /* Drawing the polygons */
  2203.  
  2204. /* HBB 980303: new functions to handle Cross back-buffer (saves
  2205.  * gp_alloc()'s) */
  2206.  
  2207. static GP_INLINE void init_Cross_store()
  2208. {
  2209.     assert(!Cross_store && !last_Cross_store);
  2210.     last_Cross_store = CROSS_STORE_INCREASE;
  2211.     free_Cross_store = 0;
  2212.     Cross_store = (struct Cross *)
  2213.     gp_alloc((unsigned long) last_Cross_store * sizeof(struct Cross),
  2214.          "hidden cross store");
  2215. }
  2216.  
  2217. static GP_INLINE struct Cross *get_Cross_from_store()
  2218. {
  2219.     while (last_Cross_store <= free_Cross_store) {
  2220.     last_Cross_store += CROSS_STORE_INCREASE;
  2221.     Cross_store = (struct Cross *)
  2222.         gp_realloc(Cross_store,
  2223.          (unsigned long) last_Cross_store * sizeof(struct Cross),
  2224.                "hidden cross store");
  2225.     }
  2226.     return Cross_store + (free_Cross_store++);
  2227. }
  2228.  
  2229. static GP_INLINE TBOOLEAN
  2230.  hl_buffer_set(xv, yv)
  2231. int xv, yv;
  2232. {
  2233.     struct Cross GPHUGE *c;
  2234.     /*HBB 961110: lclint wanted this: */
  2235.     assert(hl_buffer != 0);
  2236.     for (c = hl_buffer[xv]; c != NULL; c = c->next)
  2237.     if (c->a <= yv && c->b >= yv) {
  2238.         return TRUE;
  2239.     }
  2240.     return FALSE;
  2241. }
  2242.  
  2243. /* HBB 961201: new var's, to speed up free()ing hl_buffer later */
  2244. static int hl_buff_xmin, hl_buff_xmax;
  2245.  
  2246. /* HBB 961124: new routine. All this occured twice in
  2247.  * draw_clip_line_buffer */
  2248. /* Store a line crossing the x interval around xv between y = ya and
  2249.  * y = yb in the hl_buffer */
  2250. static GP_INLINE void update_hl_buffer_column(xv, ya, yb)
  2251. int xv, ya, yb;
  2252. {
  2253.     struct Cross GPHUGE *GPHUGE * cross, GPHUGE * cross2;
  2254.  
  2255.     /* First, ensure that ya <= yb */
  2256.     if (ya > yb) {
  2257.     int y_temp = yb;
  2258.     yb = ya;
  2259.     ya = y_temp;
  2260.     }
  2261.     /* loop over all previous crossings at this x-value */
  2262.     for (cross = hl_buffer + xv; 1; cross = &(*cross)->next) {
  2263.     if (*cross == NULL) {
  2264.         /* first or new highest crossing at this x-value */
  2265.  
  2266.         /* HBB 980303: new method to allocate Cross structures */
  2267.         *cross = get_Cross_from_store();
  2268.         (*cross)->a = ya;
  2269.         (*cross)->b = yb;
  2270.         (*cross)->next = NULL;
  2271.         /* HBB 961201: keep track of x-range of hl_buffer, to
  2272.          * speedup free()ing it */
  2273.         if (xv < hl_buff_xmin)
  2274.         hl_buff_xmin = xv;
  2275.         if (xv > hl_buff_xmax)
  2276.         hl_buff_xmax = xv;
  2277.         break;
  2278.     }
  2279.     if (yb < (*cross)->a - 1) {
  2280.         /* crossing below 'cross', create new entry before 'cross' */
  2281.         cross2 = *cross;
  2282.         /* HBB 980303: new method to allocate Cross structures */
  2283.         *cross = get_Cross_from_store();
  2284.         (*cross)->a = ya;
  2285.         (*cross)->b = yb;
  2286.         (*cross)->next = cross2;
  2287.         break;
  2288.     } else if (ya <= (*cross)->b + 1) {
  2289.         /* crossing overlaps or covers 'cross' */
  2290.         if (ya < (*cross)->a)
  2291.         (*cross)->a = ya;
  2292.         if (yb > (*cross)->b) {
  2293.         if ((*cross)->next && (*cross)->next->a <= yb) {
  2294.             /* crossing spans all the way up to 'cross->next' so
  2295.              * unite them */
  2296.             cross2 = (*cross)->next;
  2297.             (*cross)->b = cross2->b;
  2298.             (*cross)->next = cross2->next;
  2299.         } else
  2300.             (*cross)->b = yb;
  2301.         }
  2302.         break;
  2303.     }
  2304.     }
  2305.     return;
  2306. }
  2307.  
  2308.  
  2309. static void draw_clip_line_buffer(x1, y1, x2, y2)
  2310. int x1, y1, x2, y2;
  2311.      /* Draw a line in the hl_buffer */
  2312. {
  2313.     register int xv, yv, errx, erry, err;
  2314.     register int dy, nstep;
  2315.     register int ya;
  2316.     int i;
  2317.  
  2318.     if (!clip_line(&x1, &y1, &x2, &y2)) {
  2319.     /*printf("dcl_buffer: clipped away!\n"); */
  2320.     return;
  2321.     }
  2322.     if (x1 > x2) {
  2323.     errx = XREDUCE(x1) - XREDUCE(x2);
  2324.     erry = YREDUCE(y1) - YREDUCE(y2);
  2325.     xv = XREDUCE(x2) - XREDUCE(xleft);
  2326.     yv = YREDUCE(y2) - YREDUCE(ybot);
  2327.     } else {
  2328.     errx = XREDUCE(x2) - XREDUCE(x1);
  2329.     erry = YREDUCE(y2) - YREDUCE(y1);
  2330.     xv = XREDUCE(x1) - XREDUCE(xleft);
  2331.     yv = YREDUCE(y1) - YREDUCE(ybot);
  2332.     }
  2333.     if (erry > 0)
  2334.     dy = 1;
  2335.     else {
  2336.     dy = -1;
  2337.     erry = -erry;
  2338.     }
  2339.     nstep = errx + erry;
  2340.     err = -errx + erry;
  2341.     errx <<= 1;
  2342.     erry <<= 1;
  2343.     ya = yv;
  2344.  
  2345.     for (i = 0; i < nstep; i++) {
  2346.     if (err < 0) {
  2347.         update_hl_buffer_column(xv, ya, yv);
  2348.         xv++;
  2349.         ya = yv;
  2350.         err += erry;
  2351.     } else {
  2352.         yv += dy;
  2353.         err -= errx;
  2354.     }
  2355.     }
  2356.     (void) update_hl_buffer_column(xv, ya, yv);
  2357.     return;
  2358. }
  2359.  
  2360.  
  2361. /* HBB 961124: improve similarity of this routine with
  2362.  * draw_clip_line_buffer ()
  2363.  *
  2364.  * HBB 961124: changed from 'virtual' to 'do_draw', for clarity (this
  2365.  * also affects the routines calling this one, so beware!
  2366.  *
  2367.  * HBB 961124: introduced checking code, to search for the 'holes in
  2368.  * the surface' bug
  2369.  *
  2370.  * Draw a line into the bitmap (**cpnt), and optionally also to the
  2371.  * terminal
  2372.  */
  2373. static void draw_clip_line_update(x1, y1, x2, y2, do_draw)
  2374. int x1, y1, x2, y2;
  2375. TBOOLEAN do_draw;
  2376. {
  2377.     /* HBB 961110: made 'flag' a boolean variable, which needs some
  2378.      *  other changes below as well
  2379.      *
  2380.      * HBB 970508: renamed 'flag' to 'is_drawn' (reverting its meaning).
  2381.      * Should be clearer now
  2382.      */
  2383.     TBOOLEAN is_drawn;
  2384.     register struct termentry *t = term;
  2385.     register int xv, yv, errx, erry, err;
  2386.     register int xvr, yvr;
  2387.     int xve, yve;
  2388.     register int dy, nstep = 0, dyr;
  2389.     int i;
  2390.     if (!clip_line(&x1, &y1, &x2, &y2)) {
  2391.     /*printf("dcl_update: clipped away!\n"); */
  2392.     return;
  2393.     }
  2394.     if (x1 > x2) {
  2395.     xvr = x2;
  2396.     yvr = y2;
  2397.     xve = x1;
  2398.     yve = y1;
  2399.     } else {
  2400.     xvr = x1;
  2401.     yvr = y1;
  2402.     xve = x2;
  2403.     yve = y2;
  2404.     }
  2405.     errx = XREDUCE(xve) - XREDUCE(xvr);
  2406.     erry = YREDUCE(yve) - YREDUCE(yvr);
  2407.     if (erry > 0)
  2408.     dy = 1;
  2409.     else {
  2410.     dy = -1;
  2411.     erry = -erry;
  2412.     }
  2413.     dyr = dy * yfact;
  2414.     nstep = errx + erry;
  2415.     err = -errx + erry;
  2416.     errx <<= 1;
  2417.     erry <<= 1;
  2418.     xv = XREDUCE(xvr) - XREDUCE(xleft);
  2419.     yv = YREDUCE(yvr) - YREDUCE(ybot);
  2420.     (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2421.  
  2422.     is_drawn = (IS_UNSET(xv, yv) && (do_draw || hl_buffer_set(xv, yv)));
  2423.  
  2424.     /* HBB 961110: lclint wanted these: */
  2425.     assert(ymin_hl != 0);
  2426.     assert(ymax_hl != 0);
  2427.     /* Check first point in hl-buffer */
  2428.     if (xv < xmin_hl)
  2429.     xmin_hl = xv;
  2430.     if (xv > xmax_hl)
  2431.     xmax_hl = xv;
  2432.     if (yv < ymin_hl[xv])
  2433.     ymin_hl[xv] = yv;
  2434.     if (yv > ymax_hl[xv])
  2435.     ymax_hl[xv] = yv;
  2436.     for (i = 0; i < nstep; i++) {
  2437.     /* 970722: new variables */
  2438.     unsigned int xvr_temp = xvr, yvr_temp = yvr;
  2439.     if (err < 0) {
  2440.         xv++;
  2441.         xvr += xfact;
  2442.         err += erry;
  2443.     } else {
  2444.         yv += dy;
  2445.         yvr += dyr;
  2446.         err -= errx;
  2447.     }
  2448.     if (IS_UNSET(xv, yv) && (do_draw || hl_buffer_set(xv, yv))) {
  2449.         if (!is_drawn) {
  2450.         /* 970722: If this line was meant to be drawn, originally,
  2451.          * and it was just the bitmap that stopped it, then draw
  2452.          * it a bit longer (i.e.: start one pixel earlier
  2453.          */
  2454.         if (do_draw)
  2455.             (*t->move) ((unsigned int) xvr_temp, (unsigned int) yvr_temp);
  2456.         else
  2457.             (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2458.         is_drawn = TRUE;
  2459.         }
  2460.     } else {
  2461.         if (is_drawn) {
  2462.         /* 970722: If this line is only drawn because hl_buffer_set()
  2463.          * was true, then don't extend to the new position (where it
  2464.          * isn't true any more). Draw the line one pixel shorter,
  2465.          * instead:
  2466.          */
  2467.         if (do_draw)
  2468.             (*t->vector) ((unsigned int) xvr, (unsigned int) yvr);
  2469.         else {
  2470.             (*t->vector) ((unsigned int) xvr_temp, (unsigned int) yvr_temp);
  2471.             (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2472.         }
  2473.         is_drawn = FALSE;
  2474.         }
  2475.     }
  2476.     if (xv < xmin_hl)
  2477.         xmin_hl = xv;
  2478.     if (xv > xmax_hl)
  2479.         xmax_hl = xv;
  2480.     if (yv < ymin_hl[xv])
  2481.         ymin_hl[xv] = yv;
  2482.     if (yv > ymax_hl[xv])
  2483.         ymax_hl[xv] = yv;
  2484.     }
  2485.     if (is_drawn)
  2486.     (*t->vector) ((unsigned int) xvr, (unsigned int) yvr);
  2487.     return;
  2488. }
  2489.  
  2490. #define DRAW_VERTEX(v, x, y) do { \
  2491.    if ((v)->style >= 0 &&                                            \
  2492.        !clip_point(x,y)  &&                                        \
  2493.        IS_UNSET(XREDUCE(x)-XREDUCE(xleft),YREDUCE(y)-YREDUCE(ybot))) \
  2494.      (*t->point)(x,y, v->style);                                   \
  2495.    (v)->style = -1;                                                  \
  2496. } while (0)
  2497.  
  2498. /* Two routines to emulate move/vector sequence using line drawing routine. */
  2499. static unsigned int move_pos_x, move_pos_y;
  2500.  
  2501. void clip_move(x, y)
  2502. unsigned int x, y;
  2503. {
  2504.     move_pos_x = x;
  2505.     move_pos_y = y;
  2506. }
  2507.  
  2508. void clip_vector(x, y)
  2509. unsigned int x, y;
  2510. {
  2511.     draw_clip_line(move_pos_x, move_pos_y, x, y);
  2512.     move_pos_x = x;
  2513.     move_pos_y = y;
  2514. }
  2515.  
  2516. static GP_INLINE void clip_vector_h(x, y)
  2517. int x, y;
  2518.      /* Draw a line on terminal and update bitmap */
  2519. {
  2520.     draw_clip_line_update(move_pos_x, move_pos_y, x, y, TRUE);
  2521.     move_pos_x = x;
  2522.     move_pos_y = y;
  2523. }
  2524.  
  2525.  
  2526. static GP_INLINE void clip_vector_virtual(x, y)
  2527. int x, y;
  2528.      /* update bitmap, do not really draw the line */
  2529. {
  2530.     draw_clip_line_update(move_pos_x, move_pos_y, x, y, FALSE);
  2531.     move_pos_x = x;
  2532.     move_pos_y = y;
  2533. }
  2534.  
  2535. static GP_INLINE void clip_vector_buffer(x, y)
  2536.      /* draw a line in the hl_buffer */
  2537. int x, y;
  2538. {
  2539.     draw_clip_line_buffer(move_pos_x, move_pos_y, x, y);
  2540.     move_pos_x = x;
  2541.     move_pos_y = y;
  2542. }
  2543.  
  2544. static void draw(p)
  2545. struct Polygon GPHUGE *p;
  2546. {
  2547.     struct Vertex GPHUGE *v;
  2548.     struct Polygon GPHUGE *q;
  2549.     long Q;
  2550.     int i;
  2551.     int x, y;            /* point in terminal coordinates */
  2552.     register struct termentry *t = term;
  2553.     t_pnt mask1, mask2;
  2554.     long int indx1, indx2, k;
  2555.     tp_pnt cpnt;
  2556.  
  2557.     xmin_hl = HL_EXTENT_X_MAX;    /* HBB 961124: parametrized this value */
  2558.     xmax_hl = 0;
  2559.  
  2560.     /* HBB 961201: store initial values for range of hl_buffer: */
  2561.     hl_buff_xmin = XREDUCE(xright) - XREDUCE(xleft);
  2562.     hl_buff_xmax = 0;
  2563.  
  2564.     /* Write the lines of all polygons with the same id as p in the hl_buffer */
  2565.     /* HBB 961201: use next_frag pointers to loop over fragments: */
  2566.     for (q = p; (Q = q->next_frag) >= 0 && (q = plist + Q) != p;) {
  2567.     if (q->id != p->id)
  2568.         continue;
  2569.  
  2570.     /* draw the lines of q into hl_buffer: */
  2571.     v = vlist + q->vertex[0];
  2572.     TERMCOORD(v, x, y);
  2573.     clip_move(x, y);
  2574.     for (i = 1; i < q->n; i++) {
  2575.         v = vlist + q->vertex[i];
  2576.         TERMCOORD(v, x, y);
  2577.         if (q->line & (1L << (i - 1)))
  2578.         clip_vector_buffer(x, y);
  2579.         else
  2580.         clip_move(x, y);
  2581.     }
  2582.     if (q->line & (1L << (i - 1))) {
  2583.         v = vlist + q->vertex[0];
  2584.         TERMCOORD(v, x, y);
  2585.         clip_vector_buffer(x, y);
  2586.     }
  2587.     }
  2588.  
  2589.  
  2590.     /* first, set the line and point styles: */
  2591.     {
  2592.     struct lp_style_type lp;
  2593.  
  2594.     lp = *(p->lp);
  2595.     lp.l_type = p->style;
  2596.     term_apply_lp_properties(&(lp));
  2597.     }
  2598.  
  2599.     /*print_polygon (p, "drawing"); */
  2600.  
  2601.     /* draw the lines of p */
  2602.     v = vlist + p->vertex[0];
  2603.     TERMCOORD(v, x, y);
  2604.     clip_move(x, y);
  2605.     DRAW_VERTEX(v, x, y);
  2606.     for (i = 1; i < p->n; i++) {
  2607.     v = vlist + p->vertex[i];
  2608.     TERMCOORD(v, x, y);
  2609.     if (p->line & (1L << (i - 1)))
  2610.         clip_vector_h(x, y);
  2611.     else
  2612.         clip_vector_virtual(x, y);
  2613.     DRAW_VERTEX(v, x, y);
  2614.     }
  2615.     TERMCOORD(vlist + p->vertex[0], x, y);
  2616.     if (p->line & (1L << (i - 1)))
  2617.     clip_vector_h(x, y);
  2618.     else
  2619.     clip_vector_virtual(x, y);
  2620.  
  2621.  
  2622.     /* reset the hl_buffer */
  2623.     /*HBB 961110: lclint wanted this: */
  2624.     assert(hl_buffer);
  2625.     for (i = hl_buff_xmin; i <= hl_buff_xmax; i++) {
  2626.     /* HBB 980303: one part was removed here. It isn't needed any more,
  2627.      * with the global store for Cross structs. */
  2628.     hl_buffer[i] = NULL;
  2629.     }
  2630.     /* HBB 980303: instead, set back the free pointer of the Cross store: */
  2631.     free_Cross_store = 0;
  2632.  
  2633.     /* now mark the area as being filled in the bitmap. */
  2634.     /* HBB 971115: Yes, I do know that xmin_hl is unsigned, and the
  2635.      * first test is thus useless. But I'd like to keep it anyway ;-) */
  2636.     if (xmin_hl < 0 || xmax_hl > XREDUCE(xright) - XREDUCE(xleft))
  2637.     graph_error("Logic error #3 in hidden line");
  2638.     /* HBB 961110: lclint wanted these: */
  2639.     assert(ymin_hl != 0);
  2640.     assert(ymax_hl != 0);
  2641.     assert(pnt != 0);
  2642.     if (xmin_hl < xmax_hl)
  2643.     for (i = xmin_hl; i <= xmax_hl; i++) {
  2644.         if (ymin_hl[i] == HL_EXTENT_Y_MAX)
  2645.         graph_error("Logic error #2 in hidden line");
  2646.         if (pnt[i] == 0) {
  2647.         pnt[i] = (t_pnt *) gp_alloc((unsigned long) y_malloc, "hidden ymalloc");
  2648.         memset(pnt[i], 0, (size_t) y_malloc);
  2649.         }
  2650.         if (ymin_hl[i] < 0 || ymax_hl[i] > YREDUCE(ytop) - YREDUCE(ybot))
  2651.         graph_error("Logic error #1 in hidden line");
  2652.         /* HBB 970508: this shift was wordsize dependent
  2653.          * I think I fixed that
  2654.          */
  2655.         indx1 = ymin_hl[i] / PNT_BITS;
  2656.         indx2 = ymax_hl[i] / PNT_BITS;
  2657.         mask1 = PNT_MAX << (((unsigned) ymin_hl[i]) % PNT_BITS);
  2658.         mask2 = PNT_MAX >> ((~((unsigned) ymax_hl[i])) % PNT_BITS);
  2659.         cpnt = pnt[i] + indx1;
  2660.         if (indx1 == indx2)
  2661.         *cpnt |= (mask1 & mask2);
  2662.         else {
  2663.         *cpnt++ |= mask1;
  2664.         for (k = indx1 + 1; k < indx2; k++)
  2665.             *cpnt++ = PNT_MAX;
  2666.         *cpnt |= mask2;
  2667.         }
  2668.         ymin_hl[i] = HL_EXTENT_Y_MAX;
  2669.         ymax_hl[i] = 0;
  2670.     }
  2671. }
  2672.  
  2673.  
  2674. void plot3d_hidden(plots, pcount)
  2675. struct surface_points *plots;
  2676. int pcount;
  2677. {
  2678.     long Last, This;
  2679.     long i;
  2680.  
  2681. #if defined(DOS16) || defined(WIN16)
  2682.     /* HBB 980309: Ensure that Polygon Structs exactly fit a 64K segment. The
  2683.      * problem this prevents is that even with 'huge pointers', a single
  2684.      * struct may not cross a 64K boundary: it will wrap around to the
  2685.      * beginning of that 64K segment. Someone ought to complain about
  2686.      * this nuisance, somewhere... :-( */
  2687.     assert(0 == (0x1 << 16) % (sizeof(struct Polygon)));
  2688. #endif
  2689.  
  2690.     /* Initialize global variables */
  2691.     y_malloc = (2 + (YREDUCE(ytop) >> 4) - (YREDUCE(ybot) >> 4)) * sizeof(t_pnt);
  2692.     /* ymin_hl, ymax_hl: */
  2693.     i = sizeof(t_hl_extent_y) * (XREDUCE(xright) - XREDUCE(xleft) + 1);
  2694.     ymin_hl = (t_hl_extent_y *) gp_alloc((unsigned long) i, "hidden ymin_hl");
  2695.     ymax_hl = (t_hl_extent_y *) gp_alloc((unsigned long) i, "hidden ymax_hl");
  2696.     for (i = (XREDUCE(xright) - XREDUCE(xleft)); i >= 0; i--) {
  2697.     ymin_hl[i] = HL_EXTENT_Y_MAX;
  2698.     ymax_hl[i] = 0;
  2699.     }
  2700.     /* hl_buffer: */
  2701.     /* HBB 980303 new: initialize the global store for Cross structs: */
  2702.     init_Cross_store();
  2703.     i = XREDUCE(xright) - XREDUCE(xleft) + 1;
  2704.     hl_buffer =
  2705.     (struct Cross GPHUGE * GPHUGE *) gp_alloc((unsigned long) (i * sizeof(struct Cross GPHUGE *)),
  2706.                           "hidden hl_buffer");
  2707.     while (--i >= 0)
  2708.     hl_buffer[i] = (struct Cross *) 0;
  2709.  
  2710.     init_polygons(plots, pcount);
  2711.  
  2712.     /* HBB 970326: if no polygons survived the cutting away of undefined
  2713.      * and/or outranged vertices, bail out with a warning message.
  2714.      * Without this, sort_by_zmax gets a SegFault */
  2715.     if (!pfree) {
  2716.     /* HBB 980701: new code to deallocate all the structures allocated up
  2717.      * to this point. Should fix the core dump reported by 
  2718.      * Stefan A. Deutscher today  */
  2719.     if (ymin_hl) {
  2720.         free(ymin_hl);
  2721.         ymin_hl = 0;
  2722.     }
  2723.     if (ymax_hl) {
  2724.         free(ymax_hl);
  2725.         ymax_hl = 0;
  2726.     }
  2727.     if (hl_buffer) {
  2728.         free(hl_buffer);
  2729.         hl_buffer = 0;
  2730.     }
  2731.     if (Cross_store) {
  2732.         free(Cross_store);
  2733.         Cross_store = 0;
  2734.         last_Cross_store = 0;
  2735.     }
  2736.     graph_error("*All* polygons undefined or out of range, thus no plot.");
  2737.     }
  2738.     sort_by_zmax();
  2739.  
  2740.     /* sort the polygons and draw them in one loop */
  2741.     Last = -1L;
  2742.     This = pfirst;
  2743.     /* Search first polygon, that can be drawn */
  2744.     while (Last == -1L && This >= 0L) {
  2745.     This = pfirst;
  2746.     if (obscured(plist + This)) {
  2747.         Last = This;
  2748.         continue;
  2749.     }
  2750.     if (plist[This].tested != is_tested && !in_front(Last, This))
  2751.         continue;
  2752.     draw(plist + This);
  2753.     Last = This;
  2754.     }
  2755.     /* Draw the other polygons */
  2756.     for (This = plist[Last].next; This >= 0L; This = plist[Last].next) {
  2757.     if (obscured(plist + This)) {
  2758.         Last = This;
  2759.         continue;
  2760.     }
  2761.     if (plist[This].tested != is_tested && !in_front(Last, This))
  2762.         continue;
  2763.     draw(plist + This);
  2764.     Last = This;
  2765.     }
  2766.  
  2767.     /* Free memory */
  2768.     if (plist) {
  2769.     for (This = 0L; This < pfree; This++)
  2770.         free(plist[This].vertex);
  2771.     free(plist);
  2772.     plist = 0;
  2773.     }
  2774.     if (vlist) {
  2775.     free(vlist);
  2776.     vlist = 0;
  2777.     }
  2778.     if (ymin_hl) {
  2779.     free(ymin_hl);
  2780.     ymin_hl = 0;
  2781.     }
  2782.     if (ymax_hl) {
  2783.     free(ymax_hl);
  2784.     ymax_hl = 0;
  2785.     }
  2786.     if (hl_buffer) {
  2787.     free(hl_buffer);
  2788.     hl_buffer = 0;
  2789.     }
  2790.     if (Cross_store) {
  2791.     free(Cross_store);
  2792.     Cross_store = 0;
  2793.     last_Cross_store = 0;
  2794.     }
  2795. }
  2796.